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ABSTRACT 


Multiple  mesoscale  numerical  weather  simulations  are  conducted  to  evaluate  whether 
Ensemble  Sensitivity  Analysis  (ESA)  is  a  useful  tool  determining  the  sensitivity  of  a 
severe  downslope  windstorm  (DSWS)  in  complex  terrain  to  initial  conditions  and 
assimilated  observations.  A  96-member  ensemble  is  implemented  with  1.33  km  grid 
spacing.  Sensitive  regions  are  found  both  upstream  and  downstream,  based  on  a  new 
forecast  metric  that  indicates  the  potential  for  turbulence  and  strong  winds  reaching  the 
Earth’s  surface.  Approximating  the  effects  of  assimilating  a  perfect  observation  at  these 
sensitivity  locations,  then  executing  non-linear  ensemble  forecasts,  shows  that  the  linear 
approximations  in  ESA  are  reasonable.  We  analyze  the  roles  of  upstream  wind  and 
stability  structures,  and  leeside  conditions,  in  determining  the  strength  and  propagation  of 
winds  down  the  mountain  slope  and  onto  the  adjacent  plains.  Results  suggest  that  ESA  is 
a  viable  method  to  identify  observation  locations  to  improve  forecasts  of  fine-scale,  non¬ 
linear,  high-impact  events  such  as  DSWS.  Also,  fourteen  severe  DSWS  identified  by  the 
High  Wind  Alert  System  located  at  the  USAF  Academy,  CO  are  modeled  utilizing  a 
deterministic  WRF  configuration.  Analysis  shows  that  non-wave  breaking  events  account 
for  the  strongest  DSWS  and  propagate  further  away  from  the  mountains  than  breaking 
events  at  this  location. 
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I.  INTRODUCTION 


Downslope  windstorms  (DSWS)  in  regions  of  complex  terrain  are  well- 
documented,  but  also  notoriously  difficult  to  forecast  for  individual  locations  due  to 
dramatic  differences  in  observed  wind  intensity  that  can  occur  over  very  short  spatial  and 
temporal  scales.  During  certain  seasons  in  regions  of  complex  terrain,  high  winds  are 
often  the  most  disruptive  meteorological  variable  to  military  operations.  In  particular, 
aviation  missions  can  be  severely  hampered  by  high  winds  and  associated  waves,  rotors, 
and  turbulence.  One  such  location  is  at  the  United  States  Air  Force  Academy  (USAFA) 
just  north  of  Colorado  Springs,  Colorado.  In  response  to  the  disruptive  nature  of  high 
winds  to  military  operations  at  USAFA,  a  meso-net  of  twelve  automated  meteorological 
observing  stations  was  deployed  in  various  locations  along  the  Front  Range  of  the  Rocky 
Mountains  in  and  around  USAFA  in  2004.  This  network  was  aptly  named  the  High  Wind 
Alert  System  (HWAS). 

It  is  of  the  utmost  importance  to  the  military  to  determine  if  the  deployment  of 
similar  observation  networks  into  the  forward  battlespace  can  increase  the  forecast  skill 
of  high  resolution  Numerical  Weather  Prediction  (NWP)  models  in  complex  terrain, 
especially  in  data-denied  or  data-sparse  regions.  As  the  technology  and  delivery 
platforms  such  as  drones  or  UAV/UAS’s  for  these  automated  forward  deployed  systems 
are  developed,  an  important  cost  analysis  process  will  need  to  occur.  A  better 
understanding  of  the  usefulness  of  these  observations  for  characterizing  the  battlespace, 
and  understanding  these  observations  impact  on  numerical  forecasts  in  complex  terrain  is 
paramount. 

A  tool  that  can  potentially  assess  the  impact  of  observations  to  NWP  forecasts  is 
Ensemble  Sensitivity  Analysis  (ESA).  ESA  is  emerging  as  a  viable  alternative  to  adjoint 
sensitivity  and  its  potential  has  recently  been  demonstrated  by  several  studies  (Ancell  and 
Hakim  2007;  Torn  and  Hakim  2008,  2009;  Torn  2010).  ESA  can  potentially  identify 
forecast  sensitivity  to  initial  condition  and  parameter  changes.  Combined  with  an 
ensemble  data  assimilation  system,  ESA  can  help  evaluate  existing  and  proposed 

observation  network  designs.  Most  previous  work  focused  on  synoptic  scales,  and  much 

1 


less  is  known  about  ESA  for  mesoscale  predictions  with  grid  spacing  of  a  few  km. 
Because  ESA  relies  on  an  ensemble  to  sample  initial  condition  and  forecast  PDFs,  it  is 
subject  to  sampling  error  and  an  assumption  of  linear  statistics.  The  validity  of  the 
linearity  assumptions,  and  the  role  of  sampling  error  associated  with  potentially  weak 
correlations,  are  largely  unknown  for  smaller  spatial  and  faster  temporal  scales.  Because 
mesoscale  predictability  varies  according  to  the  phenomenon  of  interest  and  forecast 
metric,  the  performance  and  potential  of  ESA  will  be  a  function  of  both  phenomenon  and 
metric. 

For  decades  it  has  been  widely  accepted  within  the  meteorological  community 
that  the  forecast  accuracy  of  NWP  models  is  inextricably  linked  to  the  accuracy  of  the 
initial  conditions  of  the  forecast  model.  In  response  to  this  relationship,  one  area  of  focus 
in  research  to  improve  NWP  forecasts  has  emphasized  the  importance  of  robust  data 
assimilation  (DA)  systems  for  NWP  models.  The  goal  of  most  DA  systems  is  to  combine 
previous  model  forecasts,  current  observations,  and  the  error  statistics  of  both  to  produce 
the  most  accurate  analysis  of  the  current  state  of  the  atmosphere.  Although  not  always 
true,  the  governing  assumption  in  NWP  is  that  the  more  accurate  the  analysis,  the  more 
skillful  a  forecast  will  be. 

In  response  to  the  notion  that  a  more  accurate  analysis  of  the  atmosphere  is 
needed,  an  immense  amount  of  research  has  gone  into  producing  new  DA  systems  that 
can  handle  the  myriad  of  different  observation  types.  In  addition,  research  has  addressed 
the  question  of  what  type,  and  where,  observations  are  needed  in  the  first  place  to 
improve  analysis  accuracy.  From  this  notion  emerged  what  is  now  known  as  targeted,  or 
adaptive,  observing  aimed  at  improving  forecasts  by  taking  observations  in  regions  where 
analysis  errors  are  likely  to  be  large  and  to  grow  rapidly,  and  to  which  some  specified 
measure  of  forecast  error  is  sensitive  (Snyder  1996;  Lorenz  and  Emmanuel  1998; 
Berliner  et  al.  1999).  Movable  observing  platforms  such  as  unmanned  aircraft  are  one 
approach.  Targeted  observing  motivates  the  question:  If  the  atmosphere  could  be 
observed  at  any  place  and  time,  where  and  what  measurements  are  needed  to  improve  the 
forecast?  Although  there  has  been  much  progress  in  this  area  in  global  and  regional 
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models,  many  unresolved  questions  remain  about  the  impact  of  observations  on  the 
predictive  skill  at  fine  scale  resolutions  of  less  than  4  km. 

One  fundamental  question  that  has  to  be  addressed  in  targeted  observing  research 
is:  what  is  the  impact  of  a  particular  observation  on  the  analysis  and  thus  the  forecast?  If 
adding  a  new  observation  to  the  DA  scheme  either  doesn’t  change  the  analysis  or  if  the 
new  observation  changes  the  analysis  but  doesn’t  change  the  forecast  of  the  NWP  model, 
then  why  add  the  observation  in  the  first  place?  In  the  most  general  sense,  if  a  forecast  is 
highly  sensitive  to  initial  conditions,  then  the  forecast  is  inherently  less  predictable.  If  a 
forecast  is  not  affected  by  inaccurate  or  varied  initial  conditions,  then  that  forecast  is 
more  stable  and  thus  much  more  predictable. 

This  idea  of  sensitivity  to  initial  conditions  and  the  relation  of  the  stability  of  the 
system  to  the  limits  of  predictability  was  linked  to  the  ideas  of  chaos  and  explored  by 
Lorenz  (Lorenz  1963,  1969).  In  NWP,  much  progress  has  been  made  through  the 
realization  that  the  chaotic  behavior  of  the  atmosphere  requires  the  replacement  of  single 
deterministic  forecasts  by  ensembles  of  forecasts  with  differences  in  the  initial  conditions 
(Kalnay  2009).  In  addition,  for  an  ideal  ensemble  that  accurately  accounts  for  all  sources 
of  forecast  uncertainty,  the  verifying  truth  should  be  indistinguishable  from  the  members 
of  the  forecast  ensemble.  Also,  the  spread  of  such  an  ideal  forecast  ensemble  ought  to  be 
correlated  to  the  forecast  uncertainty  where  cases  with  large  (small)  ensemble  spread  are 
associated  with  large  (small)  forecast  uncertainty  (e.g.  Grimit  and  Mass  2007). 

Our  research  aims  to  further  the  exploration  of  forecast  sensitivity  and  forecast 
uncertainty  in  the  context  of  a  strongly  forced  mesoscale  phenomena  with  non-linear 
dynamics  by: 

•  Evaluating  the  potential  of  ESA  as  an  analysis  tool  for  determining  the 
sensitivity  of  a  severe  DSWS  in  complex  terrain  to  initial  conditions  and 
assimilated  observations. 

•  Using  ESA  to  determine  the  leeside  conditioning  versus  upstream  synoptic 
forecast  sensitivities  for  a  DSWS  to  the  initial  conditions  in  a  mature 
ensemble  forecasting  system. 

•  Evaluating  the  potential  of  ESA  to  determine  a  priori  where  and  when 
targeted  observations  can  add  value  to  a  forecast  by  reducing  ensemble 
spread  and  thus  increasing  the  reliability  of  the  forecasts. 
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A.  MOTIVATION 

DSWS  have  been  well  documented  and  modeled  for  decades.  However,  many  of 
these  studies  have  occurred  at  Boulder,  CO  or  in  the  Owens  Valley,  CA  in  the  lee  of  the 
Sierras.  These  locations  lend  themselves  to  test  and  model  theories  in  nature,  yet  in  an 
environment  that  are  similar  to  idealized  two-dimensional  ridges  used  to  develop  much  of 
the  theory.  The  geography  in  and  around  the  Pikes  Peak  region  is  more  three-dimensional 
relative  to  most  of  the  previous  studies,  and  it  follows  that  the  wind  storm  characteristics 
are  likely  more  three  dimensional  as  well. 

For  the  purposes  of  this  study,  a  severe  DSWS  is  defined  as  westerly  wind  gusts 
of  50  kt  or  greater  (NWS  severe  criteria)  recorded  by  at  least  two  HWAS  sensors.  Severe 
DSWS  have  been  recorded  regularly  using  the  HWAS  network  at  USAFA  since  its 
installation  in  2004.  A  careful  examination  of  these  events  reveals  downslope  storms  that 
are  produced  frequently  and  on  the  same  order  of  magnitude  as  found  in  Boulder  and  Ft 
Collins,  CO  to  the  north.  A  total  of  110  severe  DSWS  (Figure  1)  have  been  recorded  by 
HWAS  over  an  eight  year  period  from  2004  to  2012  (2009  excluded  due  to  data  outages) 
with  frequent  gusts  observed  greater  than  70  kt  (36  m  s'1). 
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Figure  1.  Total  number  of  severe  DSWS  broken  down  by  month  as  recorded  by  HWAS 

sensors  at  USAFA  from  2004  to  2012. 


While  the  initiation,  upstream  conditions,  and  synoptic  setup  for  severe 
downslope  wind  storms  has  been  thoroughly  examined  in  research  over  the  past  several 
decades,  less  work  has  been  done  on  defining  the  conditions  that  promote  the  propagation 
of  strong  winds  down  the  slope  and  away  from  the  mountains  into  the  valleys  or  plains. 
In  addition,  often  severe  downslope  winds  caused  by  mountain  waves  are  accompanied 
by  leeside  trapped  waves,  low-level  vortices  or  rotor,  and  flow  separation  (Doyle  2002; 
Hertensteim  and  Kuettner  2005;  Vosper  et  al.  2006;  Jiang  et  al.  2007).  Even  less  work 
has  been  done  on  describing  the  impact  these  leeside  features  have  on  the  propagation  of 
the  strong  downslope  winds  away  from  the  slope  and  into  the  valleys  or  plains. 

We  use  the  term  “leeside  conditioning”  in  our  study  to  describe  the  conditions 
prior  to  and  during  the  onset  of  high  winds  lee  of  the  mountain  slope  in  the  region  where 
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downslope  winds  are  observed  that  play  a  role  in  determining  the  distance  strong  winds 
will  move  down  the  slope.  Leeside  conditioning  can  consist  of  the  thermal  structure, 
stability  profiles,  surface  pressure  distributions,  wave  structures,  and  boundary  layer 
characteristics  which  lend  themselves  to  strong  wind  events  at  lower  elevations  separated 
from  the  mountain  slope.  Leeside  conditioning  is  often  a  mesoscale  feature  that  is  related 
to  or  aided  by  larger  synoptic  scales,  but  can  vary  at  fine  spatial  and  temporal  scales 
driven  by  the  mesoscale  details. 

Recent  studies  have  examined  the  sensitivity  of  mountain  waves  to  topography 
and  upstream  conditions  using  an  adjoint  (Doyle  2007)  and  ensemble  sensitivities 
(Reinecke  and  Durran  2009),  but  do  not  focus  on  leeside  conditioning  and  its  role  in  the 
spatial  extent  of  strong  downslope  winds.  Additionally,  studies  have  examined  leeside 
boundary  layer  flow  separations  (Jiang  et  al.  2007),  the  influence  of  inversions  (Vosper 
et  al.  2006;  Hertenstein  2009),  surface  friction  (Doyle  and  Durran  2002),  and  wind 
profiles  (Hertenstein  and  Kuettner  2005)  on  trapped  waves  and  rotors.  However,  almost 
all  these  studies  cite  the  need  for  further  examination  of  these  features  in  irregular  terrain 
and  with  the  aid  of  a  model  capable  of  performing  full  three-dimensional  simulations.  In 
summary,  the  full  exploitation  of  ESA  to  examine  the  sensitivity  of  a  forecast  to  initial 
conditions  and  assimilated  observations  has  yet  to  be  accomplished  for  a  DSWS,  and  the 
associated  leeside  conditioning,  in  more  three-dimensional  terrain.  Here  we  address  these 
issues. 
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II.  BACKGROUND 


A.  ENSEMBLE  SENSITIVITY  ANALYSIS 

1.  Overview 

A  method  that  has  shown  promise  in  identifying  locations  of  initial  condition 
sensitivity  is  Ensemble  Sensitivity  Analysis  (ESA).  ESA,  as  proposed  by  Hakim  and 
Torn  (2008),  computes  the  sensitivity  of  a  forecast  metric  (J )  with  respect  to  the  model 
initial  conditions  (  xa  )  by  utilizing  the  statistics  of  the  ensemble.  Hakim  and  Torn  used 
ESA  to  diagnose  relationships  between  the  minimum  central  pressure  of  a  mid-latitude 
cyclone  and  initial  condition  fields  such  as  geopotential  height  and  wind  speed. 

ESA  was  also  explored  by  Ancell  and  Hakim  (2007),  where  they  relate  and  show 
the  equivalence  of  ensemble  sensitivity  to  adjoint  sensitivity.  Furthermore,  Tom  and 
Hakim  (2008)  illustrated  how  ensemble  sensitivity  for  synoptic-scale  forecasts  can  be 
used  to  determine  climatological  sensitivity.  Their  experiment  demonstrated  that 
forecasts  of  sea-level  pressure  and  rainfall  over  western  Washington  were  sensitive  to  the 
upstream  mass  and  temperature  fields  (Figure  2).  Their  results  suggest  that  ESA  provides 
an  attractive  alternative  to  adjoint  sensitivity  analysis  and  may  be  useful  for  observation 
thinning  and  observation  targeting. 
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Figure  2.  Results  of  the  ESA  performed  by  Hakim  and  Torn,  showing  that  the 

sensitivity  of  24  hour  precipitation  in  a  box  over  western  Washington  State  is 
most  sensitive  to  changes  in  SLP  in  the  darkest  regions  west  of  the  box.  After 
Hakim  and  Torn  2008. Derivation  from  Linear  Regression 

A  scalar  approximation  to  the  full  ESA  can  be  derived.  Given  a  bivariate  dataset 
of  observations,  linear  regression  finds  the  straight  line  that  explains  the  linear 
relationship  of  the  sample.  This  line  is  given  by  the  equation  y  =  mx  +  b  where  y  is  the 
dependent  variable,  x  is  the  independent  variable,  m  is  the  slope  of  the  line,  and  b  is  the  y- 
intercept.  This  line  minimizes  the  sum  of  squared  differences  between  observed  values 
(the  y  values)  and  predicted  values  given  by  the  line.  Linear  regression  is  most 
appropriate  when  both  the  x  and  y  values  are  normally  distributed,  meaning  the  marginal 
distributions  are  symmetric  and  have  a  single  mode.  The  slope  of  this  best  fit  line  (m)  can 

be  expressed  by  — .  It  can  be  shown  by  using  any  bivariate  sample  that  the  slope  of  this 
best  fit  line  can  be  represented  by: 

Ay  _  cov(x,  y) 

Ax  var(.r) 

This  expression  can  be  simplified  to: 
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Ay  _  ray 
Ay  cr 


(2.2) 


where  r  is  the  correlation  coefficient  between  x  and  y,  oy  is  the  standard  deviation  of  y, 
and  ox  is  the  standard  deviation  of  x. 


ESA  in  its  most  basic  scalar  form  performs  linear  regression  to  find  the  best  fit 
line  between  two  samples  where  the  sensitivity  is  the  slope  of  this  line.  For  an  ensemble 
size  M,  the  sensitivity  of  a  scalar  forecast  metric  derived  from  the  ensemble  mean 
forecast  ( Je )  to  an  analysis  state  variable  (xa)  is  equivalent  to  Eq.  (1)  in  Tom  and  Hakim 
(2008): 

SJe  _  cov(  J  ,  x) 

Sxa  var(x) 

ESA  uses  samples  of  the  state  at  the  initial  and  forecast  time  to  estimate 
statistically  (probabilistically)  how  changes  in  the  initial  conditions  affect  the  forecast 
metric  J.  Here  J  and  x  are  1  x  M  ensemble  estimates  of  the  state  variable  (independent 
variable)  at  an  analysis  grid  point  and  forecast  metric  (dependent  variable)  at  a  defined 

location  and  time.  ESA  is  shown  to  approximate  the  adjoint  (^)  given  that  as  the  limit 
of  8Je  and  8xa  approaches  zero: 


"  '  A  (2.4) 

8xa  dx  dxa 

Switching  from  a  scalar  form  to  matrix  form,  ESA  can  be  performed  at  every  grid 
point  and  every  state  variable  represented  by  xa  that  is  a  N  x  M  matrix  where  N  is  the 
number  of  state  variables  and  M  is  the  number  of  ensemble  members. 


To  determine  the  impact  of  hypothetical  observations,  Ancell  and  Hakim  (2007) 
propose  in  their  Eq.  (26)  to  use  the  analysis  increment  and  an  ensemble  Kalman  filter 
(EnKF)  to  assimilate  observations  l..i,  and  can  be  written  as: 

Sxa  =  -x“  =  K,.  (y,°+1  -  )  (2.5) 

where 

K  =  Pa  Ht  (HPaHT  +  R)  1  (2.6) 
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and  Pf  is  the  n  x  n  ensemble-estimated  analysis  error  covariance,  H  is  the  n-length  row 
vector  containing  a  linearized  forward  operator,  and  R  is  the  p  x  p  observation  error 
variance.  It  is  important  to  note  that  the  analysis  increment  from  the  EnKF  is  consistent 
with  ensemble  statistics  only  for  the  same  EnKF  and  model.  Here  the  perturbation  to  the 
analysis  8xa  is  due  to  the  assimilation  of  a  new  observation  (yf+1).  To  determine  how 
this  perturbation  relates  to  the  forecast  metric  J,  Eq.  (2.5)  is  solved  for  8xa  substituted 
into  a  matrix  form  of  Eq.  (2.1)  yielding: 

SJ={^\  K,W„, ,-»")•  (2-7) 

\OX  J 

The  implications  of  Eq.  (2.7)  are  significant  as  it  allows,  within  the  sampling  and 
assumptions  of  ESA,  an  a  priori  assessment  of  how  a  hypothetical  observation  would 
affect  an  ensemble  mean  and  ensemble  spread  of  a  forecast  given  a  particular  model,  and 
data  assimilation  system.  Within  the  framework  presented,  ESA  is  able  to  not  only 
predict  regions  of  sensitivity  to  initial  conditions  but  also  allows  for  the  evaluation  of 
potential  targeted  observations  or  network  design  without  taking  or  assimilating  real 
observations.  This  capability  is  of  the  utmost  value  for  military  contingency  operations 
and  planning  in  data  denied  or  data  sparse  regions  where  accurate  forecasts  and 
assessment  of  forecast  reliability  are  strongly  desired. 

B.  DOWNSLOPE  WINDSTORMS 

1.  Overview 

DSWS  have  received  much  attention  by  the  meteorological  community  over  the 
past  several  decades.  Downslope  winds  are  produced  when  stably  stratified  air  encounters 
a  mountain  barrier  that  generates  gravity  waves  to  produce  strong  winds  under  certain 
conditions  in  the  lee  of  the  mountain  barrier.  Generally,  upstream  conditions  favorable 
for  strong  downslope  winds  include  the  existence  of  a  stable  layer  near  or  above 
mountain  top  and  a  less  stable  layer  above.  Winds  can  be  high  at  the  level  of  the  stable 
layer,  but  strong  winds  in  the  upper-troposphere  are  not  necessary  (Brinkman  1974). 
DSWS  have  been  given  many  names  in  different  locations  across  the  globe  and  are 
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generally  classified  into  two  types  of  winds;  foehn  winds  or  warming  winds  and  bora 
winds  or  cooling  winds. 

The  Rocky  Mountain  region  with  its  chinook  winds  is  perhaps  the  most  studied 
and  well  known  location  for  these  events.  A  chinook  is  a  foehn  wind  phenomenon  along 
the  eastern  slopes  of  the  Rocky  Mountains.  These  winds  are  much  more  frequent  and 
stronger  in  the  winter  than  in  the  summer.  The  onset  of  a  chinook  is  easily  recognized  by 
a  sharp  increase  in  temperature,  a  drop  in  relative  humidity,  the  development  of  strong 
westerly  surface  winds,  the  existence  of  a  dry-adiabatic  lapse  rate,  and  often  the 
appearance  of  wave  cloud  features  (Oard  1993).  For  a  long  time,  severe  downslope  winds 
along  the  eastern  slopes  of  the  Rocky  Mountains  were  only  attributed  to  chinook  type 
events.  However,  it  has  been  documented  that  katabatic  flows  similar  to  the  bora  of 
Europe  have  been  responsible  for  many  of  the  wind  events  along  the  Colorado  Front 
Range  (Brinkman  1974).  Some  of  Boulder’s  windstorms  have  caused  temperature  rises  of 
15-20°C  while  others  have  resulted  in  temperature  drops  of  up  to  15°C.  Despite  the 
difference  in  the  evolution  of  surface  temperature  lee  of  the  slope,  there  does  not  appear 
to  be  any  significant  dynamical  distinction  between  the  processes  responsible  for  the 
development  of  high  winds  in  cold  or  warm  events  (Durran  2003).  Regardless  of  their 
classification,  these  downslope  wind  events  can  be  intense,  damaging  to  property,  forests, 
vegetation,  military  and  aviation  operations,  and  human  life.  Wind  intensities  can  vary 
upon  small-scale  terrain  effects,  exposure,  the  synoptic  weather  patterns  and  the 
atmospheric  dynamics  present. 

Several  theories  explaining  the  dynamic  forcing  of  these  storms  have  been 
published  throughout  the  years,  with  some  disagreement.  Scorer  and  Klieforth  (1959), 
and  Aaneson  (1965)  proposed  that  the  generation  of  strong  downslope  winds  were 
primarily  attributed  to  rather  short  (generally  <  20  km  wavelength)  quasi-periodic  lee 
waves  that  propagate  downstream  of  a  mountain  range  under  conditions  that  allow 
trapping  of  wave  energy  to  create  strong  windstorms.  These  waves  form  under  conditions 
most  notably  when  the  wind  speed  increases  rapidly  with  height.  According  to  Klemp 
and  Lilly  (1975),  this  theory  proved  insufficient  to  describe  the  dynamics  of  a  DSWS  as 
observations  of  DSWS  around  Boulder,  Colorado  indicated  that  the  horizontal 

11 


wavelength  of  waves  associated  with  downslope  winds  storms  are  on  the  order  of  50-100 
km  and  vertically  propagate  freely,  leading  to  the  conclusion  that  trapped  resonant  waves 
do  not  play  an  important  role  in  this  phenomenon. 


Elias  sen  and  Palm  (1960)  demonstrated  that  when  an  upward  propagating  gravity 
wave  encounters  a  region  in  which  the  Scorer  parameter  changes  rapidly,  some  of  this 
wave  energy  can  be  directed  back  into  a  downward  propagating  wave.  The  Scorer 
parameter  is  defined  as: 


l=—r 


N  1  d2U 


(2.8) 


U2  U  dz2 

where  U( z)  is  the  flow  of  the  background  state  and  N(z)  is  the  Brunt-Vaisala  or  buoyancy 
frequency.  Extending  the  work  of  Eliassen  and  Palm  (1960),  Klemp  and  Lilly  (1975) 
suggest  that  strong  downslope  winds  occur  when  the  atmosphere  has  a  multilayer 
structure  that  produces  an  optimal  superposition  of  upward  and  downward-propagating 
waves.  They  found  that  this  most  often  occurs  when  the  tropopause  is  located  one-half 
vertical  wavelength  above  the  ground. 


Clark  and  Peltier  (1977)  further  explored  this  idea  using  a  numerical  model  and 
suggested  that  large -amplitude  waves  and  downslope  winds  are  produced  after  a 
vertically  propagating  wave  breaks;  the  energy  in  the  upward  propagating  wave  is 
trapped  beneath  its  own  level  of  “supercritical  steepening,”  producing  a  substantial 
increase  in  the  wave  amplitude.  A  series  of  studies  ensued  (Clark  and  Peltier  1984; 
Peltier  and  Clark  1979,  1984;  Clark  and  Farley  1984)  that  suggested  that  after  increasing 
in  amplitude  and  breaking,  the  wave  breaking  region  was  characterized  by  strong  mixing 
and  a  wind  reversal  in  the  cross-mountain  flow.  This  wave  induced  critical  level  provides 
a  top  to  a  resonant  layer  that  reflects  upward  propagating  waves  back  toward  the 
mountain  leading  to  very  strong  surface  winds.  Additional  numerical  studies  have 
repeatedly  verified  the  importance  of  wave  breaking  in  the  development  of  strong  surface 
winds  (Durran  1990). 


The  dynamics  governing  the  development  of  strong  downslope  winds  in  the 
atmosphere  are  similar  to  those  governing  the  rapid  increase  in  speed  that  occurs  when 

water  flowing  over  a  rock  in  a  river  undergoes  a  transition  from  a  relatively  slow  velocity 

12 


upstream  to  a  thin  layer  of  high-velocity  fluid  over  the  downstream  face.  In  such 
circumstances,  a  turbulent  hydraulic  jump  often  develops  downstream  of  the  rock  at  the 
point  where  the  high-speed  flow  decelerates  back  to  the  ambient  velocity  of  the  river 
(Durran  2003).  This  notion  has  led  researchers  to  describe  the  DSWS  on  the  basis  of 
hydraulic  jump  theory,  assuming  the  atmosphere  can  be  modeled  as  two  or  more 
neutrally  stratified  layers  with  sharp  inversions  at  each  layer  interface  (Long  1953a; 
Kuettner  1959;  Houghton  and  Isaacson  1968).  Additional  work  was  done  in  this  area 
proposing  that  strong  winds  will  occur  along  the  lee  slope  when  a  fluid  undergoes  a 
transition  from  subcritical  flow  upstream  to  supercritical  flow  over  the  mountain  (Smith 
1977;  Durran  1986;  Durran  and  Klemp  1987).  However,  it  has  been  suggested  that  the 
hydraulic  analog  is  best  applied  to  the  atmosphere  in  a  qualitative,  rather  than  a 
quantitative,  manner  (Durran  2003),  as  quantitative  application  is  made  difficult  by 
defining  a  dynamically  meaningful  Froude  number  in  vertically  unbounded  continuously 
stratified  flow  as  hydraulic  theory  bounds  the  flow  by  a  free  surface  and  hydrostatic 
balance  (Durran  1990;  Reinecke  and  Durran,  2008). 

2.  Leeside  Conditioning 

While  upstream  synoptic  forcing  is  necessary  for  the  creation  of  mountain  waves 
and  DSWS,  often  the  local  conditions  in  the  lee  of  the  mountain  barrier  can  determine  the 
severity,  impact,  and  propagation  of  the  downslope  wind  storm  down  the  slope  to  the 
adjacent  valleys  or  plains.  We  use  the  term  “leeside  conditioning”  in  our  study  to 
describe  the  conditions  prior  to  and  during  the  onset  of  high  winds  lee  of  the  mountain 
slope  in  the  region  where  downslope  winds  are  observed,  and  which  may  play  a  role  in 
determining  the  distance  strong  winds  will  move  down  the  slope.  Leeside  conditioning 
can  consist  of  the  thermal  structure,  stability  profiles,  surface  pressure  distributions,  wave 
structures,  and  boundary  layer  characteristics  which  lend  themselves  to  strong  wind 
events  at  lower  elevations  separated  from  the  mountain  slope.  Leeside  conditioning  is 
often  a  mesoscale  feature  that  varies  at  very  fine  spatial  and  temporal  scales,  but  can  be 
influenced  or  aided  by  larger  synoptic  scales. 
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a. 


Downslope  Windstorm  Leeside  Conditioning 


While  several  studies  have  considered  leeside  conditioning  during  severe 
DSWS,  the  topic  is  less  studied  than  the  upstream  conditions  which  are  thought  to  be 
important  for  the  development  or  initiation  of  DSWS.  Brinkman  (1974)  analyzed  20 
Boulder  windstorms  and  assessed  the  role  of  mean  deviations  from  normal  station 
pressure  lee  of  the  mountain  slope.  During  chinook  (warm)  events  he  found  that  lee-side 
pressure  drops  twice  as  rapidly  as  upstream  stations;  but  in  cold  or  bora  cases  pressures 
rise  twice  as  fast  in  the  upstream  stations  than  in  the  lee.  Brinkman  concludes  that  a 
strong  local  pressure  minimum  is  observed  beneath  a  warm  region  caused  by  adiabatic 
compression.  Because  this  pressure  minimum  is  caused  by  the  windstorm,  it  would  have 
little  predictive  value  and  the  extent  of  eastward  surface  winds  would  be  limited.  In 
addition,  Brinkman  noted  that  leeside  soundings  showed  a  lower  stable  layer  (inversion) 
than  upstream  soundings,  and  also  have  a  steeper  lapse  rate  for  Boulder  windstorms. 

The  depth  and  strength  of  surface-based  inversions  and  cold  pools  has 
been  identified  as  a  primary  factor  in  the  strength  of  surface  winds  during  a  downslope 
event.  Pielke  (1985)  suggested  that  in  the  absence  of  surface  heating,  cold  air  lee  of  a 
mountain  barrier  can  be  flushed  out  by  internal  gravity  waves  that  reach  sufficient 
amplitude,  and  also  by  shear-induced  turbulence  at  the  top  of  the  cold  pool.  Pielke  notes 
that  in  nature,  both  processes  are  likely  to  occur  simultaneously.  In  order  to  test  these 
theories,  Lee  et  al.  (1989)  investigated  the  influence  of  cold  pools  of  air  on  downslope 
winds  and  the  flushing  of  the  cold  air  east  of  the  mountain  barrier  using  a  two- 
dimensional  model  and  idealized  terrain.  Their  results  showed  that  cold  pools  can 
determine  whether  downslope  winds  make  it  to  the  mountain  base.  They  conclude  that 
the  top  of  the  cold  pool  acts  very  much  like  terrain  with  similar  shape  preventing  the 
penetration  of  strong  winds  to  the  mountain  base.  In  this  case,  the  core  of  maximum 
horizontal  winds  is  found  at  the  top  of  the  mountain  with  less  intense  winds  along  the  lee 
slope. 
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Jiang  and  Doyle  (2008)  conducted  a  series  of  idealized  simulations  of 
modest  downslope  events  in  the  Owens  Valley,  CA  during  the  TREX  experiment  and 
looked  at  the  effects  of  diurnal  surface  heating  on  leeside  conditioning.  They  concluded 
that  surface  heating  after  sunrise  weakens  the  downslope  winds  by  destratifying  the 
downslope  flow,  and  eventually  leads  to  the  decoupling  of  the  valley  air  from  the 
westerlies  above  the  mountaintop  by  local  noon.  Their  simulations  showed  that  a  shallow 
penetration  depth  down  the  slope  into  the  valley  exists  in  the  absence  of  surface  heating. 
Furthermore,  they  stated  that  surface  heating  in  the  lee  plays  a  dual  role;  it  tends  to 
weaken  the  downslope  flow  by  reducing  its  stratification  and  enhance  it  by  creating  a 
leeside  surface  low. 

b.  Leeside  Conditioning  Due  to  Rotor  and  Flow  Separation 

Rotor  is  defined  as  the  circulation  of  flow  about  a  horizontal  or  nearly 
horizontal  axis  that  is  usually  associated  with  flow  over  the  lee  side  of  a  barrier,  such  as  a 
mountain  range  (Glickman  2000).  Rotors  are  often  found  in  the  lee  of  a  mountain  ridge  in 
conjunction  with  high  amplitude  mountain  waves  (Holmboe  and  Kleforth  1957;  Kuettner 
1959;  Lester  and  Fingerhut  1974).  Because  a  rotor  can  exist  in  the  lee  of  mountains 
without  the  presence  of  severe  downslope  winds,  leeside  trapped  waves,  rotor,  and  flow 
separation  are  often  treated  in  the  literature  as  separate  phenomena  from  DSWS. 
However,  the  two  phenomena  are  linked  by  their  fundamental  forcing,  stratified  flow 
over  a  terrain  barrier.  This  delineation  in  most  of  the  literature  has  led  to  very  little 
research  or  linkage  between  the  convergence  zone  between  the  downslope  winds  and  the 
low-level  rotor  region. 

Kuettner  (1959)  was  probably  the  first  to  thoroughly  document  and 
research  leeside  rotors  after  experiencing  its  effects  during  many  glider  flights.  He 
summarized  the  structure  he  observed  and  proposed  in  a  diagram  shown  in  Figure  3. 
Lester  and  Fingerhut  (1974)  named  this  rotor  region  the  Lower  Turbulent  Zone  (LTZ) 
and  defined  it  as  a  highly  turbulent  region  of  nearly  neutral  stability  found  immediately  to 
the  lee  of  the  mountains  between  the  ground  and  an  elevated  stable  layer  in  which  the 
mountain  wave  motion  is  occurring.  Based  upon  the  observations  from  this  study  and 
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modeling  efforts  (Hertenstein  and  Kuettner  2005),  leeside  rotors  have  been  classified  into 
two  main  types;  lee-wave  type  and  jump  type.  In  lee  wave  types,  the  rotor  is  located 
below  lee  waves  and  contains  moderate  to  severe  turbulence.  Jump-type  rotors  extend 
well  above  the  upstream  inversion  and  contain  severe  turbulence.  Lester  and  Fingerhut 
(1974)  documented  the  varying  distances  the  first  lee  trough  can  set  up  from  the  ridge 
line  giving  a  distance  ranging  from  5  to  ~50  km.  While  much  further  research  has  been 
conducted  on  understanding  the  strength,  size,  internal  dynamics,  and  relation  to 
boundary  layer  flow  separation  for  rotors,  almost  no  progress  has  been  made  on 
determining  where  the  first  rotor  in  the  lee  of  the  ridge  line  sets  up. 
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Figure  3.  Schematic  diagram  of  flow  pattern  over  mountain  range. 

From  Kuettner  1959 
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Doyle  and  Durran  (2002)  linked  the  importance  of  lee  waves  in  the 
formation  of  rotors.  In  a  series  of  two-dimensional  numerical  simulations  they 
documented  that  the  existence  of  rotors  depends  on  the  amount  of  surface  friction;  and 
that  sufficiently  high  wave  amplitude  can  induce  a  pressure  gradient  that  can  lead  to 
boundary-layer  flow  separation.  Regions  of  recirculation  or  upstream  flow  underneath  the 
wave  crests  can  exist  as  shown  in  Figure  4.  Figure  5,  from  Doyle  and  Durran  (2002), 
shows  that  as  surface  friction  increases  the  distance  to  flow  separation  decreases. 
However,  no  further  discussion  was  given  on  what  factors  other  than  surface  friction 
might  play  a  role  in  the  locations  of  flow  separation  and  rotor  formation. 


Figure  4.  Streamlines  and  horizontal  wind  speed  (m  s'1)  for  simulations  using  (a)  free- 
slip  and  (b)  no-slip  lower  boundary  conditions.  Horizontal  wind  speed 
contours  less  than  zero  are  shaded  (from  Doyle  and  Durran  2002). 


Doyle  and  Durran  (2002)  also  discussed  the  relationship  between  the 
strength  of  downslope  winds  and  the  formation  of  lee  waves  and  rotors.  They  note  that  an 
increase  in  rotor  strength  is  often  accompanied  by  an  increase  in  the  downslope  winds 
along  the  lee  slope.  However,  they  conclude  that  it  is  the  lee-wave  amplitude  that  most 
closely  correlates  with  rotor  intensities,  rather  than  downslope  wind  strength.  The 
increase  in  the  downslope  wind  speed  is  much  smaller  than  the  increase  in  the  lee-wave 
induced  pressure  gradient.  Furthermore,  in  one  of  their  experiments  the  downslope  winds 
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were  so  strong  (exceeding  60  m  s'1)  that  no  rotors  or  reversed  flow  were  able  to  develop. 
They  conclude  that  in  the  absence  of  significant  lee  wave  activity,  rotors  do  not  need  to 
develop  even  when  the  winds  are  very  high,  but  large- amplitude  lee-waves  can  produce 
very  strong  rotors  even  when  the  downslope  winds  are  relatively  modest. 


z0  (m) 


Figure  5.  Distance  to  flow  separation  (solid,  non-dimensional)  and  pressure  drag 
(dashed)  as  a  function  of  surface  roughness  length 
(from  Doyle  and  Durran  2002). 


Jiang  et  al.  (2007)  investigated  the  relationship  between  winds  in  the  lee  of 
topography,  trapped  waves,  and  boundary-layer  separation.  They  defined  boundary-  layer 
separation  as  the  separation  of  the  boundary  layer  from  the  lee  slope  under  the  influence 
of  an  adverse  pressure  gradient.  They  described  further  that  at  the  separation  point,  the 
lowest  streamline  is  lifted  of  the  surface  with  a  recirculation  region  located  between  the 
separation  points.  They  connect  this  study  with  previous  studies,  noting  that  the 
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recirculation  region  is  also  known  as  a  rotor.  Using  a  two-dimensional  numerical  model, 
Jiang  et  al.  studied  the  effects  of  stability,  trapped  waves,  wave  breaking,  and  elevated 
inversions  on  boundary  flow  separation  lee  of  topography. 

Contrary  to  previous  studies,  such  as  Doyle  and  Durran  (2002)  and  Vosper 
(2004),  Jiang  et  al.  conclude  that  while  bottom  friction  plays  a  role  in  rotor  formation, 
surface  wind  speed  down  the  slope  of  the  terrain  is  the  most  relevant  boundary-layer 
property  that  determines  boundary-layer  separation.  Their  study  suggests  that  in  the 
presence  of  trapped  waves,  boundary-layer  separation  occurs  where  the  wave-induced 
adverse  pressure  gradient  and  the  negative  horizontal  wind  perturbation  reach  a 
maximum  and  that  boundary  layers  separation  forced  by  trapped  waves  is  governed  by  a 
normalized  ratio  of  the  vertical-velocity  maximum  to  the  surface  wind  speed.  Put  in  a 
different  way,  the  stronger  the  surface  wind  speed,  the  stronger  the  wave  forcing  needs  to 
be  to  induce  boundary  layer  separation.  While  this  study  further  enhanced  the 
understanding  of  the  conditions  during  which  boundary  layer  separation  occurs,  no 
discussion  was  given  on  what  determines  the  distance  from  the  upstream  terrain  that  the 
first  flow  separation  point  occurs. 

C.  CLIMATOLOGY  -  SEVERE  DOWNSLOPE  WINDSTORMS  AT  USAFA 

A  study  of  the  climatology  of  severe  DSWS  at  USAFA  is  conducted  analyzing 
the  HWAS  dataset  from  2004  through  2012  (2009  excluded  due  to  data  outages).  A  study 
of  this  type  is  useful  to  determine  how  the  frequency,  intensity,  and  seasonal  variation  of 
the  storms  at  USAFA  compare  to  other  commonly  studied  locations.  A  severe  DSWS  is 
defined  as  westerly  wind  gusts  of  50  kt  or  greater  (NWS  severe  criteria)  recorded  by  at 
least  two  HWAS  sensors  for  the  purposes  of  this  study.  A  total  of  110  severe  events  are 
identified  with  an  average  of  17  severe  events  per  year  (Figure  1).  Severe  low-level 
winds  also  occur  several  times  per  year  due  to  thunderstorm  outflow  or  strong  pressure 
gradients  with  winds  parallel  to  the  mountains  (winds  not  out  of  the  west)  and  these 
events  are  excluded  for  the  purposes  of  this  study.  The  complexity  of  the  terrain 
surrounding  USAFA  is  described  in  detail  in  Appendix  A.  Cross  sections  from  the 
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“Truth”  simulation  (described  in  Chapter  III)  are  given  in  Appendix  B  for  the  entirety  of 
the  30  December  2008  event  used  as  a  case  study  and  described  in  section  D. 

Severe  DSWS  at  USAFA  are  most  common  during  the  winter  months.  Peak 
frequency  of  occurrence  is  in  November  and  December,  with  a  second  peak  in  the 
months  of  March  and  April  (Figure  6).  The  second  peak  in  March  and  April  coincides 
with  the  frequent  strong  leeside  cyclogenesis  that  occurs  in  this  region  during  the  early 
spring  months,  which  may  play  a  role  in  the  development  of  leeside  troughing  that  often 
accompanies  severe  events  by  increasing  the  cross  mountain  pressure  gradient  in  this 
region.  As  shown  in  Figures  1  and  6,  severe  DSWS  occur  at  almost  any  time  of  the  year, 
with  the  exception  of  July  and  August  where  a  severe  event  has  yet  to  be  recorded. 

The  winter  of  2011/12  was  particularly  active  for  severe  DSWS  at  USAFA. 
HWAS  recorded  26  between  the  months  of  November  to  May,  with  24  of  the  26  events 
recorded  between  the  months  of  November  through  March.  An  average  of  only  6.3  days 
existed  between  severe  events  during  this  period.  The  peak  during  this  season  was  in 
January  with  six  events. 

Of  the  110  severe  DSWS  identified,  only  eight  of  these  events  were  recorded  on 
ten  or  more  of  the  twelve  sensors  and  only  two  events  were  recorded  by  all  twelve 
sensors  including  the  December  30,  2008  event  used  as  a  case  study  for  this  research.  The 
average  peak  gust  recorded  by  HWAS  for  these  eight  widespread  events  is  70  kt  (36  m  s' 
),  with  the  highest  gust  of  79  kt  (41  m  s'1)  recorded  at  the  Rampart  sensor  at  1708Z  on 
December  30,  2008.  This  event  is  discussed  in  detail  in  the  next  section  and  is  the  basis 
for  our  ESA  studies. 

HWAS  additionally  reveals  that  the  highest  winds  are  predominantly  confined  to 
the  higher  elevations  and  along  the  western-most  portions  of  USAFA  that  lie  closest  to 
the  steep  slopes  of  the  Rampart  Range.  The  Rampart  sensor  is  by  far  the  most  wind-prone 
sensor;  most  every  severe  event  was  recorded  by  this  sensor  except  for  a  few  when  this 
sensor  had  a  data  outage.  Many  more  severe  events  were  recorded  by  this  sensor  than  are 
included  in  the  110  events  cited  in  this  study  as  these  additional  events  were  only 
recorded  by  the  Rampart  sensor  and  not  by  any  of  the  other  HWAS  sensors. 
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Figure  6.  Annual  frequency  of  severe  DSWS  broken  down  by  month  as  recorded  by 

HWAS  sensors  at  USAFA  from  2004  to  2012. 


For  the  purposes  of  this  research,  the  Airfield  sensor  is  selected  to  represent  a 
valley  location  separated  from  the  slopes  of  the  Rampart  Range.  This  sensor’s  selection  is 
based  upon  its  downstream  location  almost  directly  east  of  the  higher  Rampart  sensor, 
and  because  of  the  importance  of  forecasting  winds  at  this  location  in  support  of  aviation 
activities.  An  examination  of  the  Airfield  sensor  reveals  that  only  24  of  the  110  severe 
events  (21.8%)  are  recorded  at  the  Airfield  sensor,  or  approximately  three  events  per 
year.  These  events,  although  less  frequent,  can  still  be  intense  with  gusts  as  high  as  64  kt 
(33  m  s'1).  It  is  a  primary  goal  of  this  research  to  determine  if  ESA  can  be  used  as  a  tool 
to  identify  what  initial  conditions  these  events  are  sensitive  to,  leading  to  insight  about 
the  structure  of  these  events  and  why  78%  of  severe  DSWS  observed  by  HWAS  never 
propagate  away  from  the  mountain  slope  and  move  out  into  the  adjacent  valleys  and 
plains. 
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D.  CASE  STUDY:  SEVERE  DOWNSLOPE  WINDSTORM— DECEMBER  30, 

2008 

1.  HWAS  Characterization  of  30  December  2008  Downslope  Windstorm 

The  strongest  DSWS  recorded  during  the  period  of  study  by  HWAS  was  on 
December  30,  2008.  This  event  was  recorded  by  all  twelve  HWAS  sensors  and  included 
the  strongest  gusts  ever  recorded  by  five  of  the  twelve  sensors.  This  severe  event  caused 
widespread  damage  across  Colorado  Springs  and  the  Front  Range  of  Colorado  and  is  the 
fourth  most  costly  wind  event  in  Colorado  history  with  an  estimated  $7  million  in  insured 
damage  (Rocky  Mountain  Insurance  Information  Association  2009).  Time  series  of 
winds,  wind  direction,  temperature,  dewpoint,  and  MSLP  are  shown  for  the  Rampart, 
Community  Center,  Airfield,  and  Pine  Creek  HWAS  sensors  in  Figures  7—11,  showing 
the  west  to  east  progression  during  this  severe  DSWS.  This  event  lasted  approximately 
24  hours  with  winds  ramping  up  at  the  Rampart  sensor  at  1 8Z  on  the  29th  and  continuing 
until  around  18Z  on  the  30th  (Figure  7).  At  the  Rampart  sensor,  severe  gusts  (greater  than 
50  kt)  were  observed  for  sixteen  straight  hours  from  approximately  03Z  to  19Z  on  the 
30th  with  a  max  peak  gust  of  79  kt  at  around  17Z. 
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Figure  7.  Time  series  of  winds,  wind  direction,  temperature,  dewpoint,  and  MSLP  at  the 
Rampart  sensors  during  a  severe  DSWS  on  30  December  2008. 
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Two  distinct  maxima  are  shown  in  Figure  11,  with  the  first  maximum  gust  of  62 
kt  at  05Z  (Figure  7).  After  a  relative  lull,  a  second  maximum  of  79  kt  (41  m  s'1) 
followed.  Strong  gusts  increased  at  the  Community  Center  sensor  at  the  base  of  the  slope 
at  03Z,  corresponding  to  the  beginning  of  the  first  maximum  at  the  Rampart  sensor 
(Figure  11).  Severe  gusts  were  observed  intermittently  at  the  Community  Center  sensor 
from  approximately  04Z  to  18Z  with  the  highest  gusts  of  68  kt  (35  m  s'1)  occurring  about 
1530Z  (Figure  8).  A  comparison  of  the  observed  steady  state  winds  and  “Truth”  run 
forecasts  for  the  Rampart,  Community  Center  and  Airfield  sensors  is  shown  in  Figure  1 1 . 
In  addition,  the  forecast  metric  (7)  used  in  the  ensemble  perturbation  experiments  is 
shown  for  comparison  to  the  evolution  of  severe  DSWS.  The  design  and  set-up  of  the 
WRF  “Truth”  run  forecasts,  and  development  and  calculation  of  the  forecast  metric  (7) 
are  discussed  in  Chapter  III. 


Figure  8.  Time  series  winds,  wind  direction,  temperature,  dewpoint,  and  MSLP  at  the 
Community  Center  sensor  during  a  severe  DSWS  on  30  December  2008. 
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Airfield  Sensor  30  December  2008 


Figure  9.  Time  series  winds,  wind  direction,  temperature,  dewpoint,  and  MSLP  at  the 
Airfield  sensor  during  a  severe  DSWS  on  30  December  2008. 


Pine  Creek  Sensor  30  December  2008 


Figure  10.  Time  series  winds,  wind  direction,  temperature,  dewpoint,  and  MSLP  at  the 
Pine  Creek  sensor  during  a  severe  DSWS  on  30  December  2008. 
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A  strong  cold  pool  was  present  in  the  lower  elevations  from  02Z  to  05Z  as  shown 
by  the  6°C  colder  temperatures  present  at  the  Airfield  when  compared  to  the  Community 
Center  (Figure  12).  As  winds  increased  from  the  west,  adiabatic  warming  was  observed 
first  at  the  Community  Center  sensor  from  02Z  to  04Z  while  temperatures  decreased  at 
the  Airfield  due  to  radiational  cooling  during  the  nighttime  hours.  At  04Z  the  temperature 
inversion  was  mixed  out  as  evident  from  the  rapid  warming  of  8°C  in  less  than  one  hour 
at  the  Airfield  corresponding  with  a  winds  shifting  from  light  and  variable  before  04Z  to 
westerly  and  increased  after  04Z  (Figure  9).  Winds  slowly  increased  at  the  Airfield,  but 
do  not  begin  to  prominently  increase  until  after  10Z  with  severe  winds  occurring 
intermittently  between  13Z  and  18Z  as  shown  in  Figure  9.  Severe  winds  steadily  work  to 
the  east  as  winds  increase  rapidly  from  the  west  at  1130Z  at  the  Pine  Creek  Sensor  with 
maximum  severe  gusts  of  57  kt  occurring  at  16Z  (Figure  10).  In  addition,  a  10  mb 
decrease  in  MSLP  is  observed  during  this  event  at  all  four  sensors  with  the  minimum 
pressure  occurring  between  16Z  and  17Z  corresponding  to  the  period  of  strongest  winds 
at  all  four  sensors  (Figures  7,8,9,  and  10). 
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Figure  1 1 .  Time  series  of  steady  state  wind  (kt)  at  the  Rampart,  Community  Center  and 
Airfield  sensors  during  a  severe  DSWS  on  30  December  2008  overlaid  with 
“Truth”  run  forecasts  from  WRF  simulation.  Forecast  metric  ( J)  from 
“Truth”  run  shown  with  yellow  circles  and  corresponds  to  scale  on  right  axis. 


2.  Synoptic  Conditions  of  30  December  2008  Downslope  Windstorm 

Synoptically  this  event  is  characterized  by  moderate  westerly  flow  aloft  with  the 
Grand  Junction  sounding  showing  upstream  winds  generally  between  50  -  80  kt  above 
500  mb,  and  between  20  -  30  kt  between  700  and  500  mb  (Figure  13).  This  sounding 
shows  a  relatively  high  tropopause  and  jet  for  this  time  of  the  year  at  around  200  mb.  The 
polar-  front  jet  axis  lies  to  the  north  and  dips  south  (Figure  14b)  as  part  of  an  upper-level 
trough  over  Idaho  and  northern  Utah  evident  in  the  200  mb  geopotential  height  field 
(Figure  14a)  that  extends  down  to  500  mb  (not  depicted).  A  jet  maximum  of  110  kt  lies 
along  the  Utah/Idaho  border  at  200  mb.  The  upper-level  trough  supports  a  surface  low  in 
southern  Montana  at  00Z,  which  moves  rapidly  over  the  Dakota’s  by  12Z.  This  low- 
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pressure  system,  combined  with  a  strong  high  pressure  center  over  the  four-comers 
region,  produces  a  strong  cross-mountain  surface  pressure  gradient  over  Colorado  during 
this  period  (Figure  14c  and  14d).  In  addition,  a  deep  inversion  (most  likely  due  to 
subsidence)  is  observed  in  the  Grand  Junction  sounding  from  the  surface  up  to  650  mb 
(Figure  13). 


Temp  (C) 


Rampart  Temp  Community  Center  Airfield 


Figure  12.  Time  series  of  temperature  (°C)  at  the  Rampart,  Community  Center  and 
Airfield  sensors  during  a  severe  DSWS  on  30  December  2008. 
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Figure  13.  Upstream  sounding  from  Grand  Junction,  CO  during  a  severe  DSWS  on  00Z 

30  December  2008. 
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771 


Figure  14.  a)  200  mb  geopotential  height  (00Z),  b)  200  mb  windspeed  (00Z),  c)  sea  level 
pressure  (00Z),  d)  sea  level  pressure  (12Z)  on  30  December  2008. 


Figure  15.  Water  vapor  satellite  image  from  1645Z  on  30  December  2008. 
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3.  Modeling  Study  of  the  30  December  2008  Downslope  Windstorm 

Multiple  model  runs  are  performed  for  the  30  December  2008  DSWS,  and  add 
insight  and  detail  to  the  three-dimensional  structure  and  HWAS  observations.  Two 
different  types  of  model  runs  are  used  for  this  case:  a  WRF  ensemble  run  with  96 
members  and  multiple  runs  of  a  WRF  deterministic  model.  All  models  are  run  with  the 
Advanced  Research  version  3.3.1  of  the  WRF  (ARW;  Skamarock  et  al.  2008).  Over  100 
deterministic  runs  are  performed  in  this  study,  with  the  primary  goal  of  ascertaining 
which  configuration  of  physics  packages,  number  of  vertical  levels,  time  step  needed,  and 
initialization  time  needed  to  most  accurately  simulate  this  DSWS.  Simple  verification 
metrics  for  these  runs  included  the  root  mean  square  error  with  the  HWAS  observations 
at  the  Rampart,  Community  Center,  and  Airfield  sensors.  The  information  obtained  from 
these  runs  is  primarily  used  to  determine  a  skillful  configuration  for  the  WRF  ensemble 
mns,  but  also  reveals  information  about  the  three-dimensional  structure  and  dynamics  of 
this  DSWS. 

Examination  of  the  model  runs  shows  that  a  trapped  wave  in  the  lee  of  the  highest 
terrain  of  the  Colorado  continental  divide  is  most  likely  responsible  for  the  strong 
downslope  winds  observed  at  USAFA  on  30  December  2008.  A  detailed  examination  of 
the  characteristics  of  trapped  wave  events  of  this  type  is  given  in  Chapter  V.  The  highly 
variable  and  complex  terrain  that  separates  the  Front  Range  (Rampart  Range)  and  the 
continental  divide  135  km  to  the  west  plays  a  significant  role  in  the  flow  in  this  region. 
The  downslope  flow  separates  from  the  surface  over  the  adjacent  valley,  showing  a  large- 
amplitude  rotor  and  easterly  flow  reversal  back  into  a  convergence  zone  between  the 
downslope  winds,  rotor,  and  easterly  flow  under  the  rotor.  This  is  clearly  depicted  by 
Figures  18  and  19. 
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Figure  16.  Model  terrain  height  in  meters  for  innermost  domain.  Solid  black  line 
indicates  location  of  cross  section  shown  in  Figure  17. 


Initially  at  00Z  on  30  December  2008,  the  WRF  ensemble  mean  contains 
moderate  westerly  flow  that  interacts  with  the  high  terrain  of  the  continental  divide,  but 
the  width  and  orientation  of  the  terrain  in  the  Elk  and  Sawatch  Mountains  do  not  initiate  a 
strong  wave  response.  The  Mosquito  Range  initiates  a  large  trapped  wave  that  propagates 
east  and  is  strongest  at  heights  ranging  from  4  to  6  km  above  MSL  in  height  across  the 
high  mountain  valley  of  South  Park,  and  toward  the  Rampart  Range.  This  is  shown  by 
the  higher  wind  speeds  and  tight  theta  gradient  in  this  area  shown  in  Figure  17.  The  mean 
wavelength  is  closely  linked  with  the  dominant  wavelength  of  complex  terrain  beneath 
and  culminates  with  a  downward  branch  along  the  Front  Range  adjacent  to  USAFA. 
Instead  of  the  terrain  and  low-level  inversions  reflecting  the  wave  back  upward  as  is  the 
case  to  the  west  of  the  Rampart  Range,  the  wave  dives  downward  and  accelerates  in  the 
lee  bringing  strong  winds  to  the  HWAS  sensors  and  USAFA. 
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Figure  17.  Cross  section  of  x-wind  component  (shaded)  and  potential  temperature 

(contoured)  from  west  to  east  across  Colorado  at  the  gridline  of  the  HWAS 
Rampart  sensor  for  the  30  December  2008  00Z  ensemble  mean  unperturbed 
analysis.  The  horizontal  axis  is  scaled  to  the  number  of  grid  points  from  the 
innermost  domain  with  a  resolution  of  1.33  km.  The  HWAS  sensors  and 
USAFA  are  along  the  lee  of  the  Rampart  Range  just  to  the  east  of  point  F. 
The  black  box  indicates  the  area  shown  in  the  cross  section  in  Figure  18. 
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Figure  18.  Cross  section  of  x-wind  component  (shaded)  and  potential  temperature 

(contoured)  from  west  to  east  across  the  Rampart  Range  and  adjacent  plains 
at  the  gridline  of  the  HWAS  Rampart  sensor  for  the  30  December  2008  00Z 
ensemble  mean  unperturbed  analysis.  The  horizontal  axis  is  scaled  to  the 
number  of  grid  points  from  the  innermost  domain  with  a  resolution  of  1.33 
km.  Approximate  locations  of  four  HWAS  sensors  are  marked  by  circles. 
The  red  arrow  indicates  the  horizontal  circulation  or  rotor  suggested  by  the 

cross  section. 

The  strong  downward  motion  along  the  Rampart  Range  causes  adiabatic  warming 
and  pressure  falls  in  the  lee.  The  localized  pressure  minimum  causes  a  reverse  pressure 
gradient  to  the  east  of  the  Rampart  Range,  and  easterlies  that  flow  toward  the  pressure 
minimum  are  initiated.  The  combination  of  the  upward  motion  due  to  the  trapped  wave, 
and  the  convergence  of  the  westerlies  with  the  easterlies,  creates  a  strong  upward  branch 
of  the  wave  and  boundary  layer  flow  separation  between  the  Community  Center  and 
Airfield  at  00Z  (Figure  18).  The  flow  separation  is  shown  by  the  potential  temperature 
contours  that  follow  the  surface  or  terrain  upstream  of  the  Community  Center  that 
separate  from  the  surface  and  reorient  in  the  atmosphere  at  a  level  above  the  surface.  The 
flow  separation  is  associated  with  a  large  rotor  that  extends  several  km  upward  in  the 
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Figure  19.  Cross  section  of  v-wind  component  (shaded)  and  potential  temperature 

(contoured)  from  west  to  east  across  the  Rampart  Range  and  adjacent  plains 
at  the  gridline  of  the  HWAS  Rampart  sensor  for  the  30  December  2008 
ensemble  mean  unperturbed  three  hour  forecast  valid  at  03Z.  The  horizontal 
axis  is  scaled  to  the  number  of  grid  points  from  the  innermost  domain  with  a 
resolution  of  1.33  km.  Approximate  locations  of  four  HWAS  sensors  are 

marked  by  circles. 


atmosphere  with  easterlies  underneath  the  rotor  at  the  surface  and  lowest  levels.  The  flow 
separation  reconnects  with  the  surface  to  the  east  of  Pine  Creek  where  weak  westerlies  at 
the  surface  resume.  The  1.33  km  resolution  of  this  simulation  is  not  sufficient  to  capture 
the  detailed  structures  of  the  rotor  and  flow  separation  and  only  the  gross  characteristics 
of  these  features  are  resolved  by  the  simulation  in  this  study. 

The  strong  downslope  winds  continue  to  strengthen  from  00Z  until  approximately 
12Z  and  their  eastward  extent  slowly  propagates  at  the  surface  as  shown  in  Appendix  B. 
The  convergence  zone  where  westerlies  and  easterlies  come  together,  and  boundary  layer 
separation  occurs,  is  just  west  of  the  Airfield  at  03Z  (Figure  19).  It  moves  west  of  the 
Airfield  by  06Z,  just  to  the  east  of  Pine  Creek  by  09Z,  and  just  west  of  Pine  Creek  by 

12Z.  It  can  be  clearly  seen  in  Figures  18  and  19  that  as  the  strength  of  the  downslope 
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flow  increases,  the  strength  of  the  reversed  flow  under  the  rotor  also  increases  similar  to 
results  shown  in  Jiang  et  al.  2007. 

After  12Z,  the  strongest  winds  begin  to  retreat  back  toward  the  mountains.  The 
strongest  winds  occur  just  east  of  Pine  Creek  and  over  the  Airfield  at  15Z,  and  easterlies 
replacing  westerlies  at  all  but  the  highest  terrain  are  depicted  by  18Z.  A  cold  front  passes 
through  the  area  bringing  an  end  to  the  strong  westerly  wave  activity  and  reorients  the 
flow  from  the  northeast  with  cold  air  building  at  the  lower  elevations  at  21Z  and 
deepening  to  the  top  of  the  Rampart  Range  by  00Z  on  31  December. 

The  results  from  this  ensemble  modeling  simulation  suggest  that  the  location  of 
the  boundary  layer  flow  separation  is  critical  in  determining  the  westward  extent  of  the 
severe  DSWS.  If  an  observation  location  is  east  of  the  flow  separation  and  west  of  the 
reattachment  point,  generally  weak  winds  are  found  at  the  surface  underneath  the  lee 
wave.  Winds  reverse  under  the  rotor  and  are  generally  light  and  variable  out  of  the  east. 
This  modeling  simulation  is  in  good  agreement  with  the  HWAS  observations  taken 
during  this  event.  Therefore,  identifying  the  factors  that  determine  the  distance  from  the 
upstream  ridge  that  flow  separation  occurs  is  paramount  to  determining  the  eastward 
extent  of  the  severe  winds.  A  process  for  determining  these  sensitivities  using  ESA  is 
presented  in  Chapter  III  with  results  given  in  Chapter  IV. 
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III.  DATA  AND  METHODOLOGY 


A.  OVERVIEW 

A  96-member  WRF  ensemble  system  is  implemented  for  the  30  December  2008 
DSWS  recorded  by  the  HWAS  sensors.  Details  and  characteristics  of  this  ensemble 
system  are  given  in  section  B .  Diagnosis  of  the  structure  and  detail  of  this  event  from  the 
ensemble  and  HWAS  observations  leads  to  the  development  of  a  forecast  metric  (7)  to  be 
used  in  the  ESA  calculations.  The  process  used  to  develop  the  forecast  metric  is  given  in 
section  C.  The  forecast  metric  is  tested  for  Gaussian  properties  as  described  in  section  C. 
The  ESA  calculation  is  performed  as  described  in  Chapter  II.  The  results  of  the  ESA 
calculations  are  used  to  identify  sensitive  locations  to  perturb  the  model  and  assimilate 
observations.  Methods  for  perturbing  the  ensemble  and  assimilating  observations  are 
outlined  in  section  E  and  section  F.  Results  from  these  perturbations  and  assimilated 
observations  are  given  in  Chapter  IV. 

B.  MODEL  DESIGN  AND  SETUP 

The  96-member  ensemble  system  used  in  this  research  uses  the  Advanced 
Research  version  3.3.1  of  the  WRF  (ARW;  Skamarock  et  al.  2008).  The  WRF  ensemble 
is  run  on  Cray  XE6  supercomputers  at  the  U.S.  Army  Engineer  Research  and  Development 
Center  (ERDC)  in  Vicksburg,  Mississippi.  Each  domain  in  the  model  contains  60  vertical 
rplevels  (mass)  with  the  model  top  boundary  condition  at  5  mb.  The  time  step  is  60 
seconds  for  the  outermost  domain.  The  outermost  grid  domain  (Figure  20)  uses  12  km 
horizontal  grid  spacing,  the  middle  domain  uses  4  km  horizontal  grid  spacing,  and  the 
inner  domain  uses  1.33  km  grid  spacing  and  has  a  grid  size  of  363  x  501  (north/south  x 
east/west).  The  ensemble  uses  the  parameterized  physics  configuration  shown  in  Table  1. 
Prior  to  the  ensemble  experiment,  over  100  deterministic  runs  are  analyzed  to  test  and 
determine  the  physics  suite,  number  of  vertical  levels,  and  time  step. 
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Figure  20.  WRF  nesting  configuration  used  in  this  study.  The  star  defines  the  location  of 

Colorado  Springs. 
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Table  1.  Physics  processes  used  in  the  WRF  system 


1.  Truth  Run 

To  recreate  the  conditions  that  supported  this  severe  downslope  wind  event,  an 
observing  system  simulation  experiment  (OSSE)  is  conducted.  The  experiment  begins  with 
a  single  “Truth”  WRF  5-day  simulation  run  from  26-31  December  2008.  The  WRF  truth 
simulation  is  initialized  at  0000  UTC  on  26  December  2008  using  initial  and  lateral  boundary 
conditions  provided  by  the  Climate  Forecast  System  Reanalysis  (CFSR)  (Saha  et  al.  2010). 
The  top  pressure  level  from  the  CFSR  data  is  -.2  mb.  Forecast  winds  at  three  HWAS 
sensors  from  this  truth  simulation  for  the  severe  DSWS  are  shown  in  Figure  11. 
Boundary  conditions  and  soil  conditions  are  updated  every  three  hours  with  CFSR  data. 
Cross  sections  from  this  truth  simulation  depicting  the  severe  DSWS  at  USAFA  are 
shown  in  Appendix  B. 
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2.  Ensemble  System 

Conditioning  of  the  WRF  ensemble  system  began  at  00Z  on  26  December  2008. 
Synthetic  observations  are  produced  by  interpolating  the  truth  run  to  locations  of  real 
observations  available  in  the  archived  National  Centers  for  Environmental  Prediction 
(NCEP)  PrepBUFR  data-stream.  An  ensemble  is  produced  by  assimilating  these  synthetic 
observations  with  the  Data  Assimilation  Research  Testbed  (DART).  DART,  developed 
and  maintained  at  the  National  Center  for  Atmospheric  Research  (NCAR)  is  an  open- 
source  community  resource  that  provides  well-documented  software  tools  for  data 
assimilation  education,  research  and  development  (Anderson  et  al.  2009). 

Initialization  of  the  ensemble  spread  is  completed  by  adding  spatially  consistent, 
but  random,  perturbations  drawn  from  the  default  static  WRF-Var  error  covariance. 
Assimilations  proceed  every  three  hours.  By  the  experiment  period  beginning  00Z  30 
December,  the  ensemble  initial  conditions  are  independent  from  the  random  perturbations 
and  instead  reflect  the  assimilated  synthetic  observations.  To  ensure  that  uncertainty 
continues  to  enter  the  domain,  perturbations  applied  to  lateral  boundary  conditions  are 
also  drawn  from  the  WRF-Var  static  covariance,  following  Torn  et  al.  (2006). 

An  ensemble  simulation  is  initiated  for  the  forecast  period  beginning  at  00Z  on  30 
December  running  for  24  hours  through  00Z  on  31  December  2008.  This  simulation  is 
the  control  forecast  for  all  ESA  experiments.  It  is  referred  to  as  the  no-perturbation  (NP) 
run  as  it  is  the  control  run  without  perturbations  and  is  the  base  state  of  the  system.  In 
other  experiments,  perturbations  are  applied  to  the  base  state  to  test  the  sensitivity  and 
response  of  the  ensemble  system.  The  process  for  determining  where  the  ensemble 
analysis  is  perturbed  and  details  of  how  perturbations  are  performed  are  found  in  sections 
D  and  E.  Forecasts  for  all  96  members  are  generated  every  three  hours. 

C.  FORECAST  METRIC  (J) — GAUSSIAN  TESTS 

As  described  in  Chapter  II,  the  ESA  calculation  is  a  linear  regression  between  the 
forecast  metric  (J)  and  the  initial  conditions  (xa)  of  an  ensemble.  The  underlying 
assumptions  of  linear  regression  require  all  variables  to  be  normally  distributed,  which 
implies  they  are  symmetric  and  have  a  single  mode.  The  spread  of  a  normal  distribution 
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is  given  by  its  standard  deviation  (a).  As  shown  in  Eq.  (2.2),  ESA  utilizes  the  properties 
of  normal  distributions,  their  standard  deviations,  and  correlations.  Therefore,  when  using 
ESA,  finding  forecast  metrics  (7)  that  have  normal  or  Gaussian  distributions  is  a 
requirement  for  the  linear  regression  assumptions  to  not  be  violated.  The  following 
section  describes  the  Gaussian  tests  for  normal  distributions  that  are  used  to  help  choose 
the  forecast  metric  (7). 

1.  Lilliefors  Test 

A  Lilliefors  test  is  a  modified  Kolmogorov-Smimov  test  that  determines  if  sample 
data  can  be  reasonably  fit  by  a  normal  distribution.  Lilliefors  (1967)  showed  that  the 
standard  Kolmogorov-Smirnov  test  no  longer  applies  when  certain  parameters  such  as  the 
mean  or  variance  are  not  specified  and  must  be  estimated  from  the  sample  data.  A 
Lilliefors  test  is  a  two-sided  non-parameter  goodness-of-fit  test.  This  is  in  contrast  to  the 
standard  parameter  Kolmogorov-Smimov  test,  which  is  a  one-sample  test  that  requires 
the  null  distribution  be  completely  specified.  The  Lilliefors  test  statistic  is  defined  as: 

T  =  maxx  |f(x)-5(x)|  (3.1) 

where  F(x)  is  the  hypothesized  cumulative  normal  distribution  function  with  a  mean  and 
standard  deviation  equal  to  the  mean  and  standard  deviation  of  the  sample  and  S  (x)  is  the 
sample  cumulative  distribution  function.  A  table  of  critical  values  derived  from  a  Monte 
Carlo  simulation  for  sample  sizes  less  than  1000  gives  criteria  for  goodness  of  fit.  If  the 
value  of  T  exceeds  the  critical  value  from  the  table  for  a  distribution  with  sample  size  N, 
and  a  given  level  of  significance  (a),  then  the  hypothesis  that  the  sample  is  from  a  nonnal 
distribution  is  rejected.  In  other  words,  the  number  of  observed  values  in  a  specified  bin 
of  the  hypothesized  distribution  is  far  from  what  would  be  expected  of  the  data  for  the 
hypothesized  distribution.  Hypothesis  rejection  corresponds  to  a  p-value  less  than  the 
specified  alpha  (a).  Alpha  is  set  at  .05  for  all  experiments  in  this  research,  unless 
explicitly  stated  otherwise.  The  Lilliefors  test  here  is  thus  applied  as  a  1-  a  =  .95 
confidence  interval  test. 
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2.  Choosing  a  Forecast  Metric  (J) 

A  Lilliefors  test  for  normality  was  performed  on  all  potential  candidates  for  the 
forecast  metric  J.  Several  forecast  metrics  were  initially  evaluated  at  single  locations 
before  being  rejected  in  favor  of  the  more  complex  Bulk  Richardson  metric  described  in 
the  next  section.  These  initial  locations  are  at  two  HWAS  sensor  locations;  the  Rampart 
and  Airfield  sensors.  These  locations  were  chosen  as  Rampart  is  a  high  mountain  location 
above  all  cold  pools  and  influences  of  the  valley,  flow  separation,  and  rotor  below  and  to 
the  east.  The  Airfield  sensor  is  located  below  surface-based  inversions  before  the  severe 
DSWS.  It  is  in  the  rotor  east  of  the  flow- separation  point  initially,  yet  still  experiences 
severe  criteria  winds  for  multiple  hours  on  30  December  2008.  The  exact  location  of 
these  sensors  is  shown  in  Table  2  and  Figure  51  in  Appendix  A. 

The  first  forecast  metric  tested  at  the  Rampart  sensor  was  three-dimensional  wind 
speed  (w,  v,  and  w  components).  As  shown  in  Figure  7,  the  wind  direction  changes  very 
little  at  the  Rampart  sensor  during  this  event,  allowing  a  simple  wind  speed  calculation  to 
adequately  capture  the  onset,  duration,  and  intensity  of  this  severe  DSWS  at  the  Rampart 
sensor. 

The  winds  and  temperatures  at  the  Airfield  are  far  more  varied  than  at  Rampart 
prior  to  the  DSWS,  thus  making  the  selection  of  J  much  more  complicated  at  this 
location.  The  presence  of  a  cold  pool,  inversion,  rotor,  and  flow  separation  allowed 
easterly  winds  to  persist  until  the  inversion  broke  just  after  04Z,  when  westerly  winds 
began  to  increase.  Therefore,  the  u  wind  component  is  selected  as  the  forecast  metric  at 
the  airfield  to  account  for  this  transition  from  weak  easterlies  to  strong  westerlies. 

The  three-dimensional  grid  region  over  which  the  forecast  metric  J  is  defined  is 
referred  to  as  the  “/  box.”  A  J  box  is  specified  over  each  sensor  by  the  number  of  grid 
points  east  to  west,  north  to  south,  and  vertically  to  calculate  J  on  the  finest  model 
domain.  The  function  used  to  calculate  the  forecast  metric  is  the  average  of  wind  speed 
(Rampart)  or  u  wind  component  (Airfield)  at  all  points  specified  in  a  particular  J  box. 
Various  J  box  sizes  are  tested  at  both  Rampart  and  the  Airfield  locations,  and  a  Lilliefors 
test  for  Gaussianity  is  performed  at  each  J  box  configuration.  J  box  configurations  are 
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tested  for  each  forecast  hour  (three  hour  increments)  for  the  00Z  30  December  NP 
ensemble  run  at  both  the  Rampart  and  Airfield  locations. 

While  these  forecast  metrics  are  viable  candidates  for  sensitivity  testing  due  to 
their  high  standard  deviations  and  strong  Gaussianity,  they  lack  information  about  the 
complexity  of  the  flow  characteristics  that  may  determine  the  onset  of  high  winds.  Thus, 
we  opt  to  seek  a  forecast  metric  based  on  a  more  complete  description  of  the  relevant 
physics  and  dynamic  relationships  in  the  proximity  of  the  strong  downslope  winds.  These 
additional  characteristics  include  wind  shear  and  thermal  stability. 

a.  Forecast  Metric — Modified  Bulk  Richardson  Number 

A  forecast  metric  is  developed  with  the  goal  of  capturing  the  relevant 
physics  and  structure  of  the  downslope  winds  along  the  lee  slope,  rotor  with  reversed 
surface  flow  to  the  east,  and  flow  separation  with  a  convergence  zone  between  these 
features.  A  modified  Bulk  Richardson  number  is  used  as  the  metric  to  sample  these 
relationships. 

The  Bulk  Richardson  number  is  a  dimensionless  quantity  that  uses  finite 
differences  between  levels  to  calculate  the  ratio  of  the  vertical  stability  to  the  vertical 
shear  of  the  flow.  It  is  given  by  the  equation: 

(3.2) 

0,[(A£/)2+(AV)2] 

where  g  is  the  gravitational  acceleration,  A 0V  is  the  virtual  potential  temperature 

difference  across  a  layer  of  thickness  Az ,  6V  is  the  mean  virtual  temperature  across  the 

layer,  and  A U  and  AV  are  the  changes  in  horizontal  wind  components  across  the  same 
layer  (Stull  1988). 

The  geometry  of  the  DSWS,  and  in  particular  the  angle  of  flow  separation 
from  horizontal,  leads  to  a  modified  form  of  the  Bulk  Richardson  number.  It  samples 
near-surface  buoyancy  and  wind  shear  across  the  flow  separation,  convergence  zone,  and 
strong  potential  temperature  gradient  that  existed  at  03Z  on  30  December  2008  between 
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the  Rampart  and  Airfield  locations.  This  forecast  metric  for  the  ensemble  mean  NP 
(unperturbed)  run  is  depicted  in  Figure  21  and  given  by 

J=R 

“  (3.3) 

where  9  is  the  2-m  potential  temperature,  9  is  the  2-m  mean  layer  potential  temperature, 
U  is  the  10  m  v-wind  component,  and  Az  is  the  distance  between  the  upper  and  lower 
modeled  terrain  height.  Points  1  and  2  are  determined  by  the  ensemble  mean,  with  point  2 
being  the  location  with  the  highest  local  positive  U  wind  at  10  meters  and  point  1  being 
the  location  with  the  highest  local  negative  U  wind  at  10  meters.  At  03Z,  points  1  and  2 
roughly  correspond  to  the  locations  of  the  HWAS  Rampart  sensor  and  a  location  east  of 
the  Airfield  sensor,  but  west  of  the  Pine  Creek  sensor.  A  histogram  of  the  ensemble 
members  for  the  forecast  metric  J  valid  at  03  Z  for  the  ensemble  mean  NP  (unperturbed) 
run  is  given  in  Figure  22.  A  Lilliefors  test  is  performed  on  J  at  this  valid  time  and 
confirms  that  the  metric  is  normally  distributed. 

This  modified  Bulk  Richardson  number  relates  well  to  the  strength  of  a  severe 
DSWS.  This  is  revealed  in  Figure  11,  which  shows  a  time  series  of  the  steady  state  winds 
at  the  Rampart,  Community  Center,  and  Airfield  HWAS  sensors  during  the  severe 
DSWS,  the  corresponding  forecast  from  the  truth  run  at  these  sensors,  and  the  forecast 
metric  (J)  as  calculated  from  the  truth  run.  This  plot  reveals  an  inverse  relationship 
between  the  strength  of  the  wind  event  and  the  value  of  the  Bulk  Richardson  number.  The 
Bulk  Richardson  number  steadily  decreases  in  magnitude  as  winds  increase  and  move 
from  the  higher  terrain  down  the  slope  and  to  the  adjacent  plains.  The  smallest  Bulk 
Richardson  value  (.0492)  is  calculated  at  15Z,  which  corresponds  to  the  strongest  winds 
moving  down  the  lee  slope  and  out  over  the  plains.  The  Bulk  Richardson  value  quickly 
increases  as  the  wind  event  begins  to  subside  at  18Z  and  beyond. 
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Figure  21.  Cross  section  of  x-wind  component  (shaded)  and  potential  temperature 

(contoured)  from  west  to  east  across  the  Rampart  Range  and  adjacent  plains 
at  the  gridline  of  the  HWAS  Rampart  sensor  for  the  30  December  2008 
ensemble  mean  NP  run  (unperturbed)  three  hour  forecast  valid  at  03Z.  The 
horizontal  axis  is  scaled  to  the  number  of  grid  points  from  the  innermost 
domain  with  a  resolution  of  1.33  km.  Stars  denote  location  of  points  where 
forecast  metric  J  is  calculated  using  a  modified  Bulk  Richardson  number. 
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Figure  22.  a)  Flistogram  of  the  forecast  metric  (7 )  for  96-member  NP  (unperturbed) 

WRF  ensemble  at  03Z  30  December  2008,  b)  scatterplot  of  6  vs.  (<92  -6X) ,  c) 

scatterplot  of  (U2  -f/,)2  vs.  9 ,  d)  (02  -0X)  vs.  (U:  -Ux)2 .  The  magnitude  of 
the  forecast  metric  (7)  is  shaded  in  the  scatterplots.  All  plots  are  valid  at  03Z 

30  December  2008. 

The  scatterplots  in  Figure  22  show  that  the  forecast  metric  (7)  can  be  interpreted 
in  the  traditional  Bulk  Richardson  manner.  The  highest  values  of  7  are  associated  with  a 
strong  horizontal  potential  temperature  gradient  and  weaker  wind  shear  across  the  flow 
separation  and  convergence  zone.  The  result  is  that  stability  dominates  over  wind  shear 
and  the  likelihood  of  turbulence  is  minimized.  Conversely,  low  values  of  7  are  associated 
with  strong  wind  shear  and  weak  stability;  shear  dominates  over  stability  and  the 
likelihood  of  turbulence  is  maximized. 
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R2  =  0.7362B  R2  =  0.87284 


Figure  23.  Scatterplots  of  10-m  x-wind  component  from  the  03Z  forecast  of  the  WRF 
ensemble  and  the  forecast  metric  (left)  and  2-m  potential  temperature  (right). 

In  the  context  of  the  downslope  winds,  flow  separation,  and  rotor  with  return  flow 
in  the  lower-levels,  low  values  of  J  represent  flow  situations  where  stronger  easterlies  and 
higher  potential  temperatures  are  observed  at  the  Airfield.  Conversely,  high  values  of  J 
correspond  to  weak  easterlies  at  the  Airfield  and  lower  potential  temperature.  This  is 
revealed  in  two  ways.  First,  it  is  deduced  from  the  relationships  shown  by  the  scatterplots 
in  Figure  22,  and  secondly,  from  the  fact  that  the  standard  deviation  for  both  potential 
temperature  and  the  x-wind  component  is  larger  at  the  Airfield  than  at  Rampart.  Thus,  the 
majority  of  the  variation  of  J  is  explained  by  variations  at  the  Airfield.  This  relationship 
is  also  seen  in  the  scatterplots  of  /  vs.  U2  and  /  vs.  shown  in  Figure  23,  and  is  further 
confirmed  by  the  large  coefficient  of  determination  between  these  parameters  (0.73628 
and  0.87284,  respectively). 

In  summary,  the  modified  Bulk  Richardson  number  that  is  developed  as  a  forecast 
metric  to  sample  the  relevant  physics  and  dynamic  relationships  between  the  strength  of 
the  winds  down  the  slope,  the  reverse  flow  under  the  induced  rotor,  and  flow  separation 
is  determined  to  be  Gaussian  and  behaves  in  a  manner  that  allows  for  physically 
consistent  interpretation. 
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D.  ENSEMBLE  SENSITIVITY  ANALYSIS  USING  MODIFIED  BULK 

RICHARDSON  AS  FORECAST  METRIC  (J) 

The  forecast  metric  as  described  and  calculated  in  the  previous  section  is  now 
utilized  to  perform  the  ESA  calculation  and  test  for  sensitivities  in  the  ensemble  initial 
conditions. 

1.  Addressing  Statistical  Sampling  Error 

a.  Statistical  Significance — Confidence  Test 

Given  that  the  ESA  calculation  draws  from  a  finite  ensemble,  sampling 
error  is  naturally  introduced.  In  order  to  test  for  statistical  significance  of  the  forecast 
metric  (J),  a  1  test  is  performed  at  a  95%  confidence  interval.  The  null  hypothesis  is  that 
changes  to  input  variable  xa  do  not  change  the  forecast  metric  J.  We  reject  this 
hypothesis  if  the  absolute  value  of  the  correlation  coefficient  between  xa  and  J  exceeds 
the  95%  confidence  bounds  approximated  for  the  sample  formed  by  the  96-member 
ensemble  (Tom  and  Hakim  2008). 

Probability  values  (p-values),  p[/,x°]  associated  with  each  correlation 
coefficient,  R[J,  xa]  are  computed  and  used  in  the  confidence  test.  The  p-value  is  the 
probability  of  getting  a  correlation  as  large  as  the  observed  value  by  random  chance  when 
the  true  correlation  is  zero  (Wilks  2006).  P-values  are  computed  by  transforming  the 
correlation  to  create  a  t-statistic  having  n-2  degrees  of  freedom.  These  values  are  then 
used  for  testing  the  null  hypothesis  of  no  significant  correlation.  When  using  a  95% 
confidence  interval  and  if  p[/,.V]  is  small  (less  than  0.05),  then  we  consider  the  R\J.  x"  \ 
to  be  significant  and  we  can  reject  the  null  hypothesis  (Wilks  2006).  Locations  with 
correlation  coefficients  that  fail  the  statistical  significance  test  are  rejected  and  data 
provided  from  these  points  is  withheld  from  the  sensitivity  fields. 


47 


Figure  24.  ESA  results  for  modified  Bulk  Richardson  forecast  metric  ( J)  at  model  level 
13  for  03Z  30  December  2008.  Results  for  xa=  Potential  Temperature  are  on 
the  left  and  xa  =  Specific  Humidity  are  on  the  right.  The  full  sensitivity  field 
is  shown  on  the  upper  plots  and  the  filtered  sensitivity  by  a  confidence  test  at 
a  95%  confidence  interval  is  shown  in  the  lower  plots.  Perturbation  and 
observation  assimilation  points  are  identified  within  regions  A  and  B. 

Locations  given  for  reference  are  identified  in  Figure  16. 

2.  ESA  Calculation 

The  ESA  calculation  as  given  in  eq.  2.3  is  performed  for  the  WRF  state  variables 
u,  v,  w,  T,  QVAPOR,  and  PSFC  for  each  level  and  every  grid  point  in  the  1.33  km  inner 
domain  for  three  hour  forecast  valid  at  03Z  30  December  2008.  Viable  candidates  for 
ensemble  perturbation  or  observation  assimilation  are  sought  from  the  ESA  results.  The 
sensitivity  fields  are  analyzed  level  by  level  and  grid  points  are  rank  ordered  from  the 
highest  sensitivity  values  to  the  lowest.  The  confidence  test  defined  in  the  previous 
section  is  performed  for  each  sensitivity  calculation.  Locations  with  correlation 
coefficients  that  fail  the  statistical  significance  test  and  locations  that  fail  a  Lilliefors 
Gaussianity  test  are  withheld  from  the  sensitivity  fields.  Locations  where  sensitivity 
patterns  are  found  in  multiple  variables  are  given  stronger  consideration  as  potential 
candidates  for  further  testing.  In  addition,  coherent  structures  in  both  the  horizontal  and 
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vertical  are  sought  in  the  sensitivity  fields,  as  sensitivities  that  are  linked  to  detectible 
atmospheric  features  are  desired,  not  spurious  noise  due  to  a  high  statistical  correlation. 
In  addition,  a  strong  sensitivity  upstream  and  one  downstream  of  USAFA  is  sought  to  test 
the  response  of  the  forecast  to  perturbations  and  observations  assimilated  both  upstream 
and  downstream  of  the  severe  DSWS. 

Two  regions  of  high  sensitivity  are  selected  from  the  criteria  previously  described, 
and  are  shown  in  Figure  24.  These  areas  of  sensitivity  are  found  in  several  levels  in  the 
model  but  are  largest  at  model  level  13  in  both  the  potential  temperature  and  specific 
humidity  ESA  fields.  Level  13  is  approximately  2200  m  AGL  at  both  locations  with  the 
downstream  point  at  3600  m  MSL  and  the  upstream  point  at  4150  m  MSL.  The  largest 
point  of  sensitivity  in  the  potential  temperature  in  both  areas  is  used  as  a  location  for 
introducing  perturbations  and  assimilating  observations.  The  sensitivity  point 
downstream  of  USAFA  in  region  A  is  at  location  z  =  13,  x  =  468,  and  y  =  172  at  a 
longitude  of  38.772  N  and  longitude  of  -103.09  W.  This  point  is  referred  to  as  point  1  for 
the  remainder  of  this  study.  The  sensitivity  point  upstream  of  USAFA  in  region  B  is  at 
location  z  =  13,  x  =  182,  and  y  =  162  in  the  inner  most  nested  domain  at  a  longitude  of 
38.931  N  and  longitude  of  -107.49  W.  This  point  is  referred  to  as  point  2  for  the 
remainder  of  this  study.  A  summary  of  the  location  information  and  sensitivities 
calculated  for  points  1  and  2  is  given  in  Table  2. 


z 

y 

X 

lat 

Ion 

dJ/dx 

Point  1 

13 

172 

468 

38.772 

-103.09 

-0.20031 

Point  2 

13 

162 

182 

38.931 

-107.49 

0.28948 

Table  2.  Model  grid  point  location,  latitude,  longitude  and  sensitivity  ( dJ/dx )  for 
point  1  and  point  2.  The  units  for  clJ/dx  are  per  Kelvin. 


E.  PERTURBATION  TESTS 

Perturbation  tests  are  conducted  at  each  sensitivity  location  to  examine  whether 
ESA  is  skillful  at  identifying  sensitivity  locations.  Perturbations  are  introduced  in  two 
ways.  The  first  method  introduces  perturbations  at  our  defined  sensitivity  points  and 

spreads  the  perturbation  using  a  linear  regression  of  the  ensemble  analysis  as  sampled  by 
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the  NP  ensemble  run.  The  second  method  utilizes  an  Ensemble  Kalman  Filter  (EnKF) 
and  assimilates  observations  using  full  data  assimilation  methods.  Both  methods  are 
implemented  in  order  to  compare  the  effectiveness  of  one  versus  the  other.  The 
perturbation  method  is  less  computationally  expensive  as  it  only  performs  single  variable 
linear  regression,  and  it  does  not  reduce  the  spread  of  the  analysis,  but  shifts  the 
distribution  in  the  direction  of  the  perturbation.  In  contrast,  the  full  data  assimilation 
method  is  more  similar  to  a  full  multivariate  regression  and  allows  not  only  a  shift  of  an 
analysis  distribution,  but  also  a  reduction  in  the  spread  of  the  ensemble  at  the  points  of 
observation  assimilation  in  the  analysis.  These  distinct  differences  motivate  our  decision 
to  test  both  methods. 


1.  Perturbation  Method 


Separate  ensemble  WRF  runs  are  generated  with  perturbations  introduced  first  at 
sensitivity  point  1  and  then  at  point  2  in  a  second  ensemble  run.  These  ensemble  runs  are 
referred  to  as  PI  (perturbation  at  eastern/downstream  sensitivity  point)  and  P2 
(perturbation  at  western/upstream  sensitivity  point)  for  the  remainder  of  this  experiment. 
The  perturbation  is  introduced  by  increasing  the  ensemble  mean  potential  temperature 
value  for  xa  by  a  value  taken  from  the  WRF  truth  run  at  that  location.  The  effect  of  this 
perturbation  at  a  single  point  is  spread  through  to  other  model  variables  in  the  horizontal 
and  vertical  through  linear  regression  between  the  analysis  variable  at  the  perturbation 
point  and  the  analysis  variables  from  the  original  NP  ensemble  run.  The  intent  of  the 
linear  regression  is  to  reduce  model  shock  in  the  subsequent  analysis.  Following 
Reinecke  (2008)  and  Torn  and  Hakim  (2009)  the  ensemble  mean  perturbation  to  each  ith 
analysis  state  variable  can  be  written  as: 


p  a 
X  —  X; 


CO  v(  ■<,■<) 
var(.<) 


-a 


(3.4) 


where  xf  is  the  ensemble  mean  perturbation,  x“  is  the  analysis  variable  from  the  original 


unperturbed  ensemble,  xas  is  the  analysis  variable  from  the  perturbation  point,  and  a  is 
the  scalar  perturbation  made  at  the  individual  grid  point  for  all  i-state  variables  given  by 
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a  =  x\  -  x“ 


(3.5) 


where  x*  is  the  perturbation  variable  at  the  perturbation  point  from  the  truth  run,  and  x‘‘  is 

the  perturbation  variable  at  the  perturbation  point  from  the  ensemble  mean  unperturbed 
analysis  (NP  run). 

This  perturbation  method  ensures  the  ensemble  variability  is  left  unchanged  while 
the  mean  is  shifted.  This  perturbation  is  applied  to  all  three  domains  and  then  run  forward 
to  show  a  forecast  change  due  to  the  perturbation.  The  results  of  these  perturbations  are 
shown  in  the  next  chapter. 

F.  DATA  ASSIMILATION  METHODS 

1.  Covariance  Localization 

Due  to  sampling  error  inherent  in  the  ensemble  system,  it  is  essential  to  devise  a 
method  to  limit  the  horizontal  and  vertical  magnitude  of  the  covariance  between  xa  and 
J.  Typically  a  spatial  function  p  weights  the  data  assimilation  regression  to  express 
confidence  in  the  regression.  This  spatial  function,  p  is  typically  a  compactly  supported 
function  similar  in  shape  to  a  Gaussian,  1.0  at  the  observation  location  (Gaspari  and  Cohn 
1999).  Decreasing  covariance  according  to  pis  called  covariance  localization,  and  helps 
prevent  the  creation  of  spurious  covariances  at  long  distances  arising  from  sampling 
error.  While  this  process  can  introduce  some  imbalances  to  the  state  (Houtekamer  and 
Mitchell  2001),  it  is  a  “practical  necessity  for  EnKF  NWP  applications”  (Hamill  2013). 

2.  Observation  Assimilation 

An  observation  using  the  temperature  from  the  truth  simulation  is  assimilated 
using  a  two  part  process  outlined  by  Anderson  (2003).  An  update  increment  is  computed 
for  each  prior  ensemble  estimate  of  the  observation  variable  by  applying  a  scalar 
ensemble  filter.  Second,  a  linear  regression  of  the  prior  ensemble  sample  of  each  state 
variable  on  the  observation  variable  is  performed  to  compute  update  increments  for  each 
state  variable  ensemble  member  from  corresponding  observation  variable  increments. 
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Equation  2.5  defines  how  matrix  data  assimilation  can  be  performed  with  a  single 
observation.  The  Kalman  gain  matrix  from  eq.  2.5  is  now  defined  as 

K  =  p*Pa  HT(H/?*PaHT +R)  1  (3.6) 

where  p  represents  the  covariance  localization  matrix,  Pa  the  background  error  variance 
matrix,  *  the  Hadamard  or  Schiir  product  between  the  two  arguments,  H  the  linearized 
forward  operator,  and  R  the  observation  error  covariance  matrix.  Equation  2.5  can  be 
written  in  scalar  form  for  a  single  scalar  additional  observation  as 


u  u  . 

XM=Xi  +  Px 


.2  ‘  t  .  2 

CT  +  CT 

yy  yy  o 


/  o  a\ 

(yi+ 1  -yt ) 


(3.7) 


with  p  representing  the  localization  covariance  function  between  the  arguments,  a 


xy 


the  covariance  between  the  arguments,  <J~  the  variance  between  the  arguments,  cr  the 


observation  variance,  y°+l  the  observation  value,  and  y°  the  value  from  the  background 

at  the  observation  location.  Anderson  (2003)  showed  that  introducing  an  observation  to 
the  analysis  in  this  scalar  form  is  permissible  and  equivalent  to  the  matrix  form  because 
our  DART  analysis  is  performed  in  serial.  Assimilating  a  single  observation  in  this 
manner  yields  the  same  result  as  if  it  were  included  in  the  observation  set  at  analysis  time 
in  the  original  assimilation.  This  data  assimilation  procedure  allows  for  the  contribution 
of  one  single  observation  to  an  analysis  and  thus  a  forecast  to  be  measured. 


Three  separate  WRF  ensemble  runs  are  performed  by  1)  assimilating  temperature 
from  the  WRF  truth  as  an  observation  only  at  point  1,  2)  assimilating  temperature  from 
the  WRF  truth  as  an  observation  only  at  point  2,  and  3)  assimilating  temperature  from  the 
WRF  truth  at  both  points  1  and  2  simultaneously.  These  three  ensemble  runs  are  referred 
to  as  DAI,  DA2,  and  DAJ,  respectively,  for  the  remainder  of  this  study.  The  results  of 
these  runs  are  presented  in  the  next  chapter. 
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IV.  ESA  RESULTS 


A.  OVERVIEW 

ESA  is  performed  using  the  methodology  described  in  Chapter  III.  Two 
sensitivity  regions  are  identified  (Figure  24),  with  sensitivity  point  1  coming  from  region 
A  and  sensitivity  point  2  coming  from  region  B.  Five  separate  ensemble  runs  are 
generated  (PI,  P2,  DAI,  DA2,  DAJ)  as  described  in  the  previous  chapter.  The  sizes  of  the 
perturbations  for  PI  and  P2  are  given  as  5x  in  Table  3  and  the  values  of  the  temperature 
observations  from  the  truth  run  assimilated  for  DAI,  DA2,  and  DAJ  are  given  in  Table  4. 


Theta 

Mean  NP 

Theta  Truth 

5 

Point  1 

307.5800 

307.7629 

0.1829 

Point  2 

305.9969 

306.1230 

0.1261 

Table  3.  Potential  temperature  (K)  at  points  1  and  2  from  the  ensemble  mean  of  the 
NP  run  and  from  the  truth  run  valid  at  00Z  30  December  2008.  The 
difference  between  the  truth  run  and  NP  mean  is  shown  in  the  far  right 

column  as  5x. 


Temperature 
Mean  NP 

Temperature 
Truth  Run 

5  (Theta) 

Point  1 

273.4546 

273.6215 

0.6118 

Point  2 

266.3875 

266.4537 

0.1894 

Table  4.  Temperature  (K)  at  points  1  and  2  from  the  ensemble  mean  of  the  NP  run 
and  from  the  truth  run  valid  at  00Z  30  December  2008.  The  corresponding 
perturbation  in  potential  temperature  after  assimilation  between  the  truth 
run  and  NP  mean  is  shown  in  the  far  right  column  as  8x. 
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<J> 

Expected  8/ 

Actual  8/ 

Percent 

NP 

0.2130 

0.0620 

- 

- 

- 

PI 

0.1953 

0.0639 

-0.0366 

-0.0177 

48.28% 

P2 

0.2180 

0.0586 

0.0365 

0.0050 

13.73% 

DAI 

0.1928 

0.0819 

-0.1225 

-0.0202 

16.47% 

DA2 

0.2256 

0.0609 

0.0548 

0.0126 

23.00% 

DAJ 

0.1956 

0.0766 

-0.0687 

-0.0174 

25.31% 

Table  5.  The  mean,  standard  deviation,  expected  change,  actual  change  and  percent 
of  predicted  change  of  the  forecast  metric  J  for  each  of  the  WRF  ensemble 
runs.  All  values  are  valid  at  the  forecast  time  of  03Z  30  December  2008. 


Figure  25.  Plot  of  expected  SJ  (x-axis)  vs.  actual  SJ  (y-axis)  for  each  ensemble  run  at  the 
forecast  time  of  03Z  30  December  2008.  The  red  line  is  a  perfect  1  to  1 
relationship  for  reference.  This  plot  corresponds  to  the  data  given  in  Table  5. 


ESA  identifies  two  sensitivity  regions.  While  these  regions  only  need  a  statistical 
correlation  to  show  sensitivity,  a  physical  interpretation  is  desired  to  gain  insight  into  the 
effects  perturbations  may  have  on  the  physical  structures  in  the  atmosphere.  This  also 
highlights  the  utility  of  ESA  to  detect  these  structures  from  an  independent  approach 
(statistical  correlation  only),  possibly  leading  to  the  identification  and  location  of 
important  dynamic  relationships  in  the  atmosphere  that  have  not  been  previously 
identified.  In  this  chapter,  tests  of  the  sensitivity  points  are  first  presented  corresponding 

to  each  of  the  experiments.  Then  a  physical  interpretation  is  offered  in  section  B,  a  total 
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energy  analysis  of  each  experiment  is  given  in  section  C,  and  a  method  for  accomplishing 
linearity  tests  with  fewer  ensemble  members  is  given  in  section  D. 

1.  Linearity 

A  linearity  test  is  needed  to  examine  responses  to  perturbations  propagated  by  the 
non-linear  model,  compared  to  the  response  predicted  by  the  linear  ESA.  This  tests  the 
validity  of  ESA.  Since  we  are  linearizing  about  the  ensemble  mean,  it  is  possible  for  the 
ensemble  to  evolve  in  a  non-linear  manner,  but  if  the  perturbations  are  sufficiently  small, 
then  linearity  may  hold  for  a  discreet  amount  of  time.  All  results,  unless  otherwise  noted, 
are  for  a  verification  time  for  J  of  03Z  on  30  December  2008.  As  shown  in  Table  5  and 
Figure  25,  every  experiment  moves  the  ensemble  mean  for  J  in  the  correct  direction  as 
predicted  by  the  sensitivity  ( clJ/clx )  given  in  Table  2  and  shown  in  Figure  24.  However,  as 
shown  in  Figure  25,  ESA  overestimates  the  impact  of  introducing  a  perturbation  or 
assimilating  an  observation  in  every  experiment. 

In  general,  experiments  at  point  1  (downstream)  have  a  larger  actual  SJ  and 
impact  the  forecast  more  than  experiments  at  point  2  (upstream).  This  is  partially  due  to 
the  larger  8x  introduced  for  DAI.  The  PI  ensemble  run  moves  J  as  a  percentage  of 
actual  change  to  predicted  change  more  than  any  other  perturbation.  This  PI  run  moves  J 
over  48%  toward  the  expected  SJ  as  predicted  by  the  ESA  performed  using  the  NP  run  as 
shown  in  Table  5.  The  expected  8J  is  calculated  by  multiplying  the  sensitivity  value 
(i clJ/dx )  from  a  perturbation  point  by  the  size  of  the  perturbation  (Sx)  introduced  at  the 
same  point.  While  the  perturbations  introduced  are  very  small,  even  these  small 
perturbations,  identified  by  ESA,  impact  the  forecast  metric  J,  showing  the  strong 
dependence  of  mountain  waves  on  the  initial  conditions. 

2.  Atmospheric  Structure 

Examination  of  the  perturbation  structures  imposed  by  either  linear  regression,  or 

data  assimilation  of  observations  at  analysis  time,  leads  to  physical  insight.  In  the  case  at 

hand,  ESA  successfully  identifies  the  importance  of  the  vertical  structure  of  atmospheric 

stability  and  the  vertical  wind  profile  both  upstream  and  downstream  of  the  DSWS. 

Figure  26  shows  the  vertical  profiles  of  potential  temperature,  Brunt- Vaisala  frequency, 
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and  v-wind  component  at  the  downstream  (top)  and  upstream  points  (bottom).  These 
profiles  reveal  that  the  sensitivity  region  identified  by  ESA  at  model  level  13  is  at  the  top 
of  a  stable  layer  in  the  atmosphere  at  both  points  1  and  2.  The  higher  values  of  Brunt- 
Vaisala  frequency  just  below  level  13  in  Figure  26  (b  and  f)  clearly  show  this  stable 
layer.  In  a  physical  sense,  introducing  a  perturbation  of  higher  potential  temperature,  (PI, 
P2)  or  temperature  (DAI,  DA2)  at  these  locations  has  the  effect  of  sharpening  and 
lowering  the  top  of  this  stable  layer. 

Figure  27  shows  how  the  vertical  structure  of  potential  temperature  changes  in  the 
initial  conditions  of  the  unperturbed  (NP)  run  from  west  to  east  across  the  downstream 
point.  A  strong  stable  layer,  approximately  identified  by  the  green  colors  in  the  cross 
section  is  shallowest  just  lee  of  the  front  range  and  begins  to  widen  or  deepen  moving 
toward  the  east.  A  hypothetical  increase  in  potential  temperature  at  the  top  of  this  stable 
layer  over  the  east  would  lead  to  a  flattening  of  the  potential  temperature  ridge  between 
2.5  and  3.5  km  in  the  eastern  third  of  Figure  27  over  the  eastern  plains  of  Colorado. 

Figure  28  shows  the  vertical  structure  of  potential  temperature  in  the  initial 
conditions  of  the  NP  run  from  west  to  east  across  the  upstream  point.  A  stable  layer, 
approximately  identified  by  the  green  colors  in  the  cross  section,  is  just  upstream  of  the 
central  Colorado  Rockies.  A  hypothetical  increase  in  potential  temperature  at  the  top  of 
this  stable  layer,  just  to  the  west  and  roughly  at  the  level  of  the  highest  terrain,  could 
modulate  the  frequency  of  the  wave  response  downstream.  This  could  be  accomplished 
by  decreasing  the  height  of  the  inversion  which  could  change  how  much  energy  is 
trapped  and  the  frequency  at  which  it  oscillates.  Different  dynamics  and  locations  of  the 
highest  winds  and  wave  structures  to  the  east  would  follow.  This  bolstered  stable  layer 
could  lead  to  an  increased  trapping  effect  in  this  layer  of  the  atmosphere,  possibly  leading 
to  a  stronger  trapped  response  in  the  wave  propagating  in  the  lower  atmosphere  across  the 
Rockies  or  downstream  across  the  plains. 
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(K),  b)  Brunt- Vaisala  frequency  (s'1),  c)  a;- wind  component  (m  s'1),  d)  Scorer 
parameter  and  at  point  2  (bottom)  of  e)  Potential  Temperature  (K),  f)  Brunt- 
Vaisala  frequency  (s'1),  g)  x-wind  component  (m  s'1),  h)  Scorer  parameter. 
These  plots  are  from  the  NP  ensemble  mean  valid  at  00Z  30  December  2008. 
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Figure  27.  Cross  section  from  west  to  east  (grid  points  300-500)  at  the  gridline  of  point 
1,  showing  ensemble  mean  potential  temperature  from  the  NP  run  analysis. 
Black  arrow  denotes  the  location  of  the  introduced  perturbation  (point  1). 
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Figure  28.  Cross  section  from  west  to  east  at  the  gridline  of  point  2  (grid  points  ISO- 
375),  showing  ensemble  mean  potential  temperature  from  the  NP  run 
analysis.  Black  arrow  denotes  the  location  of  the  introduced  perturbation 

(point  2). 


3.  Forecast  Certainty  and  Ensemble  Spread 

While  it  is  not  a  primary  goal  of  the  perturbation/data  assimilation  experiments 
here  to  create  a  more  accurate  forecast,  but  rather  to  test  the  sensitivity  of  the  ensemble 
system  to  changes  in  initial  conditions,  it  is  worthwhile  to  comment  on  changes  in  the 
certainty  of  the  forecast  due  to  the  perturbations/data  assimilation.  The  most  likely 
forecast  of  a  well-developed  ensemble  system  is  the  ensemble  mean  for  a  single  mode 
distribution.  The  standard  deviation  of  the  distribution  is  often  used  to  access  ensemble 
spread  and  the  certainty  of  the  forecast.  If  the  standard  deviation  decreases,  the  forecast  is 
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expected  to  be  more  certain.  If  the  standard  deviation  increases,  the  forecast  is  expected 
to  be  less  certain. 


Figure  29.  Histograms  of  perturbation  potential  temperature  at  point  lfrom  the  analysis 
for  the  a)  NP  run,  b)  PI  run,  c)  DAI  run,  and  d)  DAJ  run. 
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Figure  30.  Histograms  of  perturbation  potential  temperature  at  point  2  from  the  analysis 
for  the  a)  NP  run,  b)  P2  run,  c)  DA2  run,  and  d)  DAJ  run. 


<  0  > 

a 

NP  Run  Point  1 

7.582 

0.0789 

NP  Run  Point  2 

5.9978 

0.0606 

PI  Run  Point  1 

7.7649 

0.0789 

P2  Run  Point  2 

6.1239 

0.0606 

DAI  Run  Point  1 

8.1913 

0.0198 

DA2  Run  Point  2 

6.1863 

0.0146 

DAJ  Run  Point  1 

8.1857 

0.0198 

DAJ  Run  Point  2 

6.1791 

0.0146 

Table  6.  Table  showing  the  mean  and  standard  deviation  of  perturbation  potential 
temperature  at  points  1  and  2  from  the  analysis  for  every  ensemble  run. 


Introducing  perturbations  and  assimilating  observations  at  points  1  and  2  have 

different  effects  on  the  standard  deviation  at  these  points  in  the  analysis.  The  standard 

deviation  and  shape  of  the  distribution  does  not  change  when  introducing  a  perturbation 

(PI  and  P2  runs).  This  is  seen  when  comparing  histogram  A  to  histogram  B  in  Figures  29 

and  30.  The  standard  deviation  of  each  is  the  same,  but  the  distribution  is  shifted  in  the 
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direction  of  the  perturbation.  In  contrast,  the  standard  deviation  decreases  and  the 
distributions  shift  when  assimilating  an  observation  (DAI,  DA2,  and  DAJ  runs).  This  is 
seen  when  comparing  histograms  C  and  D  to  histogram  A  in  Figures  29  and  30.  The  data 
assimilation  process  decreases  the  standard  deviation  at  the  observation  point  and  shifts 
the  distribution  in  the  direction  of  the  observation.  This  is  caused  by  assimilating  an 
observation  with  zero  observation  error  in  our  experiments.  The  mean  and  standard 
deviation  of  perturbation  potential  temperature  at  the  analysis  time  is  shown  in  Table  6 
for  every  ensemble  run  and  corresponds  to  the  histograms  shown  in  Figures  29  and  30. 

It  is  expected  that  the  decrease  in  the  analysis  standard  deviation  will  lead  to  a 
decreased  standard  deviation  in  J  at  the  forecast  time.  However,  this  is  not  always  the 
case  for  our  experiments.  After  the  perturbations/data  assimilation,  the  03Z  forecast 
standard  deviation  of  J  increases  in  all  ensemble  runs  that  include  the  downstream  point 
(PI,  DAI,  DAJ)  and  decreases  for  all  runs  that  include  the  upstream  point  (P2,  DA2). 
While  run  PI  responds  with  a  more  accurate  change  in  J  than  any  other  run,  at  the  same 
time  its  spread  increases  as  measured  by  standard  deviation  as  shown  in  Table  5. 
Conversely,  runs  perturbed  at  point  2  respond  with  a  less  accurate  change  in  J ,  but  their 
spread  decreases  as  measured  by  standard  deviation. 

A  possible  explanation  to  why  the  spread  of  J  increases  at  the  forecast  time,  when 
it  is  expected  to  decrease,  is  that  the  mean  flow  advects  the  atmospheric  features 
containing  the  lower  spread  in  the  analysis  downstream  as  the  forecast  progresses.  This 
has  the  effect  of  moving  features  with  decreased  spread  toward  the  location  of  J  for  the 
upstream  perturbation/assimilation  point  and  moving  features  farther  away  from  J  for  the 
downstream  perturbation/assimilation  point.  This  is  supported  by  the  cross  sections 
shown  in  Figures  31,  32,  and  33. 
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Figure  31.  Cross  Sections  of  aDA2  ~<JNP  from  west  to  east  at  the  gridline  of  the  Rampart 
sensor  (grid  points  150-375)  for  a)  analysis  potential  temperature,  b)  03Z 
forecast  of  potential  temperature,  c)  analysis  Jt-wind  component,  d)  03Z 
forecast  of  x-wind  component.  The  arrow  denotes  the  longitude  of  the 
observation  insertion  (point  2),  and  the  star  denotes  the  location  of  USAFA 

and  the  forecast  metric  J. 
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Figure  32.  Cross  Sections  of  <JDAl  ~<JNP  from  west  to  east  at  the  gridline  of  the  Rampart 
sensor  (grid  points  300-500)  for  a)  analysis  potential  temperature,  b)  03Z 
forecast  of  potential  temperature,  c)  analysis  x-wind  component,  d)  03Z 
forecast  of  x-wind  component.  The  arrow  denotes  the  longitude  of  the 
observation  insertion  (point  1),  and  the  star  denotes  the  location  of  USAFA 

and  the  forecast  metric  J. 
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Figure  33.  Cross  Sections  of  aDAJ  -  crNP  from  west  to  east  at  the  gridline  of  the  Rampart 

sensor  (grid  points  300-500)  for  a)  analysis  potential  temperature,  b)  03Z 
forecast  of  potential  temperature,  c)  analysis  x-wind  component,  d)  03Z 
forecast  of  v-wind  component.  The  arrow  denotes  the  longitude  of  the 
observation  insertion  (point  1),  and  the  star  denotes  the  location  of  USAFA 

and  the  forecast  metric  J. 


Figure  31  shows  a  cross  section  of  the  difference  between  the  DA2  run  standard 
deviation  and  the  NP  run  standard  deviation  ( aDA2  -  aNP )  for  potential  temperature  (top), 

jt-wind  component  (bottom),  at  analysis  time  (left),  and  at  the  three  hour  forecast  (right). 
Figure  31a  shows  a  region  of  decreased  spread  in  potential  temperature  associated  with 
the  observation  location.  This  area  of  decreased  spread  evolves  into  a  complicated  wave 
structure,  but  clearly  moves  to  the  east  in  the  three  hour  forecast  as  shown  in  Figure  31b. 
Similarly,  decreased  spread  is  seen  in  the  x-wind  component  at  analysis  time  just  above 
the  observation  location  and  moves  to  the  east  by  the  three  hour  forecast  (Figure  31 
bottom). 
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Figure  32  shows  a  cross  section  of  the  difference  between  the  DAI  run  standard 
deviation  and  the  NP  run  standard  deviation  ( <7 ,m  -  aNP )  for  potential  temperature  (top), 

.i'-w'ind  component  (bottom),  at  analysis  time  (left),  and  at  the  three  hour  forecast  (right). 
Figure  32a  shows  a  region  of  decreased  spread  for  potential  temperature  in  the  analysis 
associated  with  the  location  of  the  assimilated  observation.  This  region  is  not  observed  in 
the  three  hour  forecast  and  based  on  the  evidence  from  the  DA2  run  of  the  decreased 
spread  features  advecting  downstream,  it  is  presumed  that  the  decreased  spread  regions  in 
DAI  advect  to  the  east  and  out  of  the  domain.  This  potentially  explains  why  the  standard 
deviation  increases  at  the  forecast  time  of  03Z  for  J  in  each  run  that  has  a 
perturbation/assimilation  point  downstream  of  the  forecast  metric’s  location.  Further 
evidence  for  this  explanation  is  presented  in  section  C  of  this  chapter  which  shows  from  a 
total  energy  perspective  that  perturbation  energy  from  the  downstream  location 
propagates  out  of  the  domain  more  quickly  than  perturbations  introduced  at  the  upstream 
location.  Similar  results  as  shown  in  the  DAI  and  DA2  runs  are  found  simultaneously  in 
the  DAJ  mn  that  is  shown  in  Figure  33. 
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Figure  34.  Cross  sections  of  jc-wind  component  (shaded)  and  potential  temperature 

(contoured)  from  west  to  east  across  the  Rampart  Range  and  adjacent  plains 
at  the  gridline  of  the  Rampart  sensor  for  the  30  December  2008  ensemble 
mean  NP  run  (top  left),  PI  run  (top  right),  and  difference  between  PI  and  NP 
(bottom)  for  the  three  hour  forecast  valid  at  03  Z. 
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Figure  35.  a)  Histogram  of  the  forecast  metric  (J )  for  PI  96-member  WRF  ensemble  at 
03Z  30  December  2008,  b)  scatterplot  of  ~6  vs.  ( 02  —0X) ,  c)  scatterplot  of 

(U2  -Ux)2  vs.  6 ,  d)  (02  -0t  )  vs.  (U2  -Ux)2 .  The  magnitude  of  the  forecast 
metric  (7)  is  shaded  in  the  scatterplots.  All  plots  are  for  03Z  30  December 

2008. 


B.  FORECAST  METRIC  IMACT 

1.  PI  Ensemble  Run 

The  PI  ensemble  run  with  the  downstream  perturbation  at  point  1,  over  the 
eastern  plains  of  Colorado,  decreases  the  forecast  metric  7  as  expected.  The  sensitivity 
prediction  ( clJ/dx )  correctly  indicates  that  an  increase  in  potential  temperature  at  point  1 
corresponds  with  a  decrease  in  7.  The  magnitude  of  the  shift  is  over  48%  of  what  ESA 
predicts  (Table  5).  As  shown  in  Figure  34,  the  decrease  in  7  is  due  to  a  slight  increase  of 
the  downslope  winds  and  a  slightly  larger  increase  in  the  easterly  flow  beneath  the  rotor 
over  the  adjacent  plains,  thus  increasing  the  wind  shear  acting  on  the  stable  layer.  In 
addition,  the  perturbation  of  potential  temperature  causes  a  slight  shift  to  the  east  of  the 
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flow  separation  and  convergence  zone,  moving  stronger  westerlies  about  2-3  km  further 
to  the  east  than  in  the  NP  run. 


Figure  36.  a)  Histogram  for  the  difference  in  the  forecast  metric  ( J )  between  PI  and  NP 
runs,  b)  scatterplot  of  the  difference  between  the  PI  and  NP  runs  for  9  vs. 

(02  -  6X) ,  c)  scatterplot  of  the  difference  between  the  PI  and  NP  runs  for 

(U2  -U j)2  vs.  6 ,  d)  scatterplot  of  the  difference  between  the  PI  and  NP  runs 
for  (<92  -6X)  vs.  (U2  -Ux)2 .  The  magnitude  of  the  difference  between  the  PI 
and  NP  runs  for  the  forecast  metric  (7)  is  shaded  in  the  scatterplots.  All  plots 
are  valid  for  03Z  30  December  2008. 

The  changes  observed  in  the  PI  run  are  physically  consistent  in  terms  of  what  is 
expected  when  J  decreases.  Since  J  is  a  modified  Bulk  Richardson  number,  a  decrease  in 
J  should  be  due  to  increased  shear  and  decreased  stability  causing  the  likelihood  of 
turbulence  to  increase.  This  can  be  seen  in  Figures  34,  35,  and  36  valid  for  the  three  hour 
forecast.  The  scatterplots  in  Figure  35  are  consistent  with  the  NP  runs  in  the 
interpretation  of  the  forecast  metric.  The  distributions  of  the  shear  and  stability  terms  in  J 
allow  for  a  consistent  physical  interpretation  of  the  results.  The  highest  values  of  J  are 
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associated  with  a  strong  horizontal  potential  temperature  gradient  and  weaker  wind  shear 
across  the  flow  separation  and  convergence  zone.  The  result  is  that  stability  dominates 
over  wind  shear  and  the  likelihood  of  turbulence  is  minimized.  Conversely,  low  values  of 
J  are  associated  with  strong  wind  shear  and  weak  stability;  shear  dominates  over  stability 
and  the  likelihood  of  turbulence  is  maximized.  Figure  36  shows  that  a  majority  of 
members  have  decreased  with  respect  to  J  due  to  both  increased  shear  and  decreased 
stability.  Some  members  have  an  increase  with  respect  to  J  due  to  decreased  shear  and 
increased  stability.  Their  numbers  are  fewer  and  their  change  is  smaller  compared  to  the 
members  where  J  decreases.  Overall  the  ensemble  mean  of  J  decreases. 

Figure  37  (top)  shows  the  perturbations  introduced  in  the  PI  run  at  point  1  for  all 
model  levels,  and  shows  the  small  potential  temperature  increase  introduced  by  adding  a 
perturbation  at  level  13.  The  perturbation  introduced  at  level  13  leads  to  greater  changes 
in  potential  temperature  in  the  upper-troposphere  and  lower  stratosphere,  as  also  shown 
by  the  cross  section  in  Figure  38.  This  results  from  a  strong  covariance  between  point  1 
and  the  upper-troposphere  in  the  ensemble  analysis  used  for  regressing  the  perturbation. 
The  perturbation  is  spread  in  the  lower  troposphere  east  to  west  along  the  potential 
temperature  ridge  (Figure  32),  with  an  effect  of  flattening  this  ridge  over  the  eastern  third 
of  the  cross  section  and  sharpening  the  vertical  transition  from  more  stable  to  less  stable 
air  in  this  region. 

Figure  37  (bottom)  reveals  an  inverse  relationship  between  potential  temperature 
and  specific  humidity  at  point  1.  The  sensitivity  analysis  shown  in  Figure  24  suggests  that 
increases  in  potential  temperature  should  have  similar  effects  as  decreasing  the  specific 
humidity  at  the  same  point.  This  relationship  holds  after  the  perturbation  is  introduced,  as 
increasing  potential  temperature  has  the  effect  of  lowering  specific  humidity  at  the 
downstream  point  (point  1). 
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Figure  37.  Perturbation  of,  a)  potential  temperature  and,  b)  specific  humidity  from  the  PI 

run  analysis  at  point  1  as  function  of  model  level  (x-axis). 
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Figure  38.  Cross  section  of  potential  temperature  difference  between  the  PI  and  NP  run 
analysis  from  west  to  east  at  the  gridline  of  point  1  (grid  points  300-500). 
Black  arrow  denotes  the  location  of  the  introduced  perturbation  (point  1). 
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Figure  39.  a)  Histogram  of  the  forecast  metric  (J )  for  P2  96-member  WRF  ensemble  at 
03Z  30  December  2008,  b)  scatterplot  of  9  vs.  ( 02  scatterplot  of 

(U2  -Ux)2  vs.  6 ,  d)  (02  ~0{)  vs.  (U2  -Ux)2 .  The  magnitude  of  the  forecast 

metric  ( J)  is  shaded  in  the  scatterplots.  All  plots  are  valid  for  03Z  30 

December  2008. 
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Figure  40.  a)  Histogram  for  the  difference  in  the  forecast  metric  ( J )  between  P2  and  NP 
runs,  b)  scatterplot  of  the  difference  between  the  P2  and  NP  runs  for  0  vs. 
(<92  -0X)  ,c)  scatterplot  of  the  difference  between  the  P2  and  NP  runs  for 

(U2  -t/j)2  vs.  9 ,  d)  scatterplot  of  the  difference  between  the  P2  and  NP  runs 

for  (02  Gy )  vs.  (U2  —  C/j ) 2 .  The  magnitude  of  the  difference  between  the  P2 

and  NP  runs  for  the  forecast  metric  (7)  is  shaded  in  the  scatterplots.  All  plots 
are  valid  for  03Z  30  December  2008. 


2.  P2  Ensemble  Run 

Changes  in  P2  due  to  the  introduction  of  a  perturbation  of  potential  temperature 
from  the  truth  run  at  the  upstream  point  (point  2)  are  almost  not  discernible  in  the  mean 
statistics.  Cross  sections  comparing  the  ensemble  mean  of  P2  to  NP  at  the  forecast  time 
and  location  of  J  are  virtually  identical  (not  shown).  However,  changes  in  the  distribution 
and  to  individual  members  are  apparent.  Figure  40  shows  that  a  majority  of  members 
exhibit  very  little  change  and  the  slight  increase  in  the  P2  J  mean  is  driven  by  relatively 
few  members  that  decrease  shear  and  increase  in  stability.  In  addition,  the  inverse 
relationship  between  potential  temperature  and  specific  humidity  predicted  by  the 
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Potential  Temp  Perturbation  (*C)  Specific  Humidity  Perturbation  (gkg*1) 


Figure  41 .  Perturbation  of  a)  potential  temperature  and  b)  specific  humidity  for  the  P2 
run  analysis  at  point  2  as  function  of  model  level  (x-axis). 

sensitivity  analysis  (Figure  41b)  holds  true  after  the  introduction  of  a  potential 
temperature  perturbation  at  point  2. 

Figure  41a  shows  the  perturbations  introduced  in  the  P2  run  in  the  column 
through  point  2,  and  shows  the  small  increase  in  potential  temperature  introduced  by 
adding  a  perturbation  at  level  13.  The  perturbation  introduced  at  level  13  leads  to  stronger 
downward  branches  of  the  downstream  trapped  waves.  This  is  seen  by  the  location  of 
increased  potential  temperature  east  of  point  2,  and  also  by  comparing  the  location  of  the 
waves  in  Figure  28  and  the  location  of  potential  temperature  increases  in  Figure  42.  The 
increased  sharpness  in  the  upstream  stable  layer  has  the  effect  of  increasing  reflection 
downward  of  the  downstream  trapped  waves  in  the  P2  run. 
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Figure  42.  Cross  section  of  ensemble  mean  potential  temperature  difference  between  the 
analysis  for  the  P2  and  NP  runs  from  west  to  east  at  the  gridline  of  point  2 
(grid  points  150-375).  The  black  arrow  denotes  the  location  of  the  introduced 

perturbation  (point  2). 

As  noted  in  the  overview  to  this  chapter,  the  spread  of  J  as  measured  by  standard 
deviation  decreases  in  the  P2  run  (Table  5).  This  is  seen  graphically  when  comparing 
histograms  between  Figures  22  and  39  for  the  NP  and  P2  runs.  The  histogram  for  the  P2 
run  shows  a  much  sharper  distribution  than  the  NP  histogram,  with  the  highest  bar  in  the 
plot  having  more  than  20  members  in  its  range.  This  suggests  increased  certainty  in  the 
forecast  of  J  at  03Z  in  the  P2  run. 
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3. 


DAI  Ensemble  Run 


The  DAI  ensemble  run  has  the  largest  change  in  J  of  any  run  due  to  the  large  5x 
introduced  through  assimilation.  However,  it  only  moves  J  slightly  more  than  16%  of 
what  ESA  predicts  it  should  move.  ESA  indicates  that  an  increased  potential  temperature 
at  the  downstream  point  causes  a  decrease  in  J.  This  is  revealed  in  the  cross  section  of  the 
ensemble  mean  in  the  region  of  J  (Figure  43)  and  in  scatterplots  (Figures  44  and  45).  The 
mean  of  J  decreases  as  the  majority  of  members  move  in  a  negative  direction.  The 
negative  SJ  results  from  increased  shear,  decreased  stability,  and  increased  mean 
potential  temperature.  The  increase  in  shear  is  primarily  due  to  a  2.5  m  s-1  increase  in  the 
easterlies  in  the  return  flow  under  the  rotor.  Additionally,  the  westerlies  and  the  region  of 
flow  separation  shift  slightly  to  the  east  causing  an  increase  of  over  5  m  s'1  in  the  region 
just  west  of  the  Airfield  between  the  DAI  and  NP  runs  in  the  ensemble  mean. 

Figure  46  (top)  shows  the  perturbations  introduced  in  the  DAI  run  at  point  1,  and 
reveals  the  potential  temperature  increase  introduced  by  assimilating  the  temperature 
from  the  truth  run  at  level  13.  Assimilating  this  temperature  observation  at  level  13  also 
leads  to  greater  potential  temperature  in  the  upper-troposphere  and  lower  stratosphere, 
(also  shown  by  the  cross  section  in  Figure  47).  The  magnitude  of  the  changes  from  data 
assimilation  is  significantly  larger  in  the  DAI  than  the  PI  case,  most  likely  due  to  a  larger 
5x  at  point  1  after  assimilation  of  the  temperature  observation.  At  most  locations  the 
perturbations  in  DAI  are  2-3  times  larger  than  the  PI  run.  The  implications  of  this  in 
terms  of  total  energy  are  explored  in  the  next  section.  Similar  to  PI,  the  perturbation  is 
spread  in  the  lower  troposphere  east  to  west  along  the  potential  temperature  ridge 
identified  in  Figure  37.  This  has  the  effect  of  flattening  this  ridge  over  the  eastern  third  of 
the  cross  section  and  sharpening  the  vertical  transition  from  more  stable  to  less  stable  air 
in  this  region.  Fastly,  Figure  46  (bottom)  reveals  an  inverse  relationship  between 
potential  temperature  and  specific  humidity  at  point  1  for  the  DAI  run,  consistent  with  all 
other  runs  and  suggested  by  the  original  ESA  prediction. 
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Figure  43.  Cross  sections  of  v-wind  component  (shaded)  and  potential  temperature 
(contoured)  from  west  to  east  across  the  Rampart  Range  and  adjacent  plains 
at  the  gridline  of  the  Rampart  sensor  for  the  30  December  2008  ensemble 
mean  NP  run  (top  left),  DAI  run  (top  right),  and  difference  between  DAI  and 
NP  (bottom)  for  the  three  hour  forecast  valid  at  03  Z.  All  plots  are  valid  at 

03Z  30  December  2008. 
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Figure  44.  a)  Histogram  of  the  forecast  metric  (7 )  for  DAI  96-member  WRF  ensemble 
at  03Z  30  December  2008,  b)  scatterplot  of  9  vs.  ( 02  -  6X) ,  c)  scatterplot  of 

(U2  U ! )2  vs.  6 ,  d)  (02  ~6{)  vs.  (U2  -f/,)2 .  The  magnitude  of  the  forecast 
metric  (7)  is  shaded  in  the  scatterplots.  All  plots  are  valid  at  03Z  30 

December  2008. 
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Figure  45.  a)  Histogram  for  the  difference  in  the  forecast  metric  (7 )  between  DAI  and 
NP  runs,  b)  scatterplot  of  the  difference  between  the  DAI  and  NP  runs  for  9 
vs.  (02  ,  c)  scatterplot  of  the  difference  between  the  DAI  and  NP  runs 

for  (U2  -  Ux)2  vs.  0 ,  d)  scatterplot  of  the  difference  between  the  DAI  and  NP 

runs  for  ( 02  -0X)  vs.  (U2  -Ux)2 .  The  magnitude  of  the  difference  between 

the  DAI  and  NP  runs  for  the  forecast  metric  (7)  is  shaded  in  the  scatterplots. 
All  plots  are  valid  at  03  Z  30  December  2008. 
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Figure  46.  Perturbation  of  a)  potential  temperature  and  b)  specific  humidity  for  the  DAI 
run  analysis  at  point  2  as  function  of  model  level  (x-axis). 
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Figure  47.  Cross  section  of  potential  temperature  difference  between  the  analysis  for  the 
DAI  and  NP  runs  from  west  to  east  at  the  gridline  of  point  1  (grid  points 
300-500).  The  black  arrow  denotes  the  location  of  the  introduced 
perturbation  (point  1). 

4.  DA2  Ensemble  Run 

The  results  of  the  DA2  run  are  very  similar  to  the  findings  in  the  P2  run,  except 
the  magnitude  of  the  changes  are  greater  when  comparing  DA2  with  P2.  The  DA2  run 
verifies  better  than  the  P2  run  in  terms  of  how  much  J  actually  changes  when  compared 
to  the  predicted  change  by  ESA.  Qualitatively  the  patterns  and  shifts  in  the  various  fields 
are  almost  identical  for  the  DA2  run  as  the  P2  run.  The  shifts  in  J  continue  to  be 
consistent  with  the  Bulk  Richardson  reasoning,  as  the  small  increase  in  J  is  caused  by 
both  increases  in  stability  and  decreases  in  shear  in  the  majority  of  the  members  (Figures 
48  and  49). 
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Figure  48.  a)  Histogram  of  the  forecast  metric  ( J )  for  DA2  96-member  WRF  ensemble 
at  03Z  30  December  2008,  b)  scatterplot  of  9  vs.  ( 02  -0X),  c)  scatterplot  of 

(U2  -t/j)2  vs.  6 ,  d)  (02  -6X)  vs.  (U2  -t/j)2 .  The  magnitude  of  the  forecast 
metric  (J)  is  shaded  in  the  scatterplots.  All  plots  are  valid  at  03Z  30 

December  2008. 

5.  DAJ  Ensemble  Run 

The  DAJ  run  has  the  same  temperature  observations  as  assimilated  at  point  1  in 
DAI  and  point  2  in  DA2,  but  jointly  assimilates  them  into  the  same  run.  Overall,  the  DAJ 
run  resembles  very  closely  the  results  obtained  from  the  DAI  run.  Cross  sections  from 
the  DAJ  run  are  not  shown  and  the  reader  is  referred  to  cross  sections  in  section  3  from 
the  DAI  run  as  the  structures  are  almost  identical.  The  perturbation  at  point  1  dominates 
the  domain  and  causes  changes  to  the  forecast  that  are  much  larger  in  magnitude  than 
changes  caused  by  the  perturbation  at  point  2,  especially  in  the  area  in  and  around  the 
forecast  metric  J  at  USAFA.  This  is  due  to  the  perturbation  associated  with  the 
downstream  condition  being  more  linearly  related  to  the  downslope  winds  three  hours 
later.  The  influence  of  observation  assimilation  at  point  2  is  seen,  but  mostly  upstream 
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over  the  higher  terrain  and  much  smaller  than  the  influences  caused  by  assimilating  the 
temperature  observation  from  the  truth  run  at  point  1. 


Figure  49.  a)  Histogram  for  the  difference  in  the  forecast  metric  (7 )  between  DA2  and 
NP  runs,  b)  scatterplot  of  the  difference  between  the  DA2  and  NP  runs  for  9 
vs.  (02  ~9\) ,  c)  scatterplot  of  the  difference  between  the  DA2  and  NP  runs 

for  (U2  -  Ux)2  vs.  0 ,  d)  scatterplot  of  the  difference  between  the  DA2  and  NP 
runs  for  {02  -6X)  vs.  (U2  -Ux)2 .  The  magnitude  of  the  difference  between 
the  DA2  and  NP  runs  for  the  forecast  metric  (7)  is  shaded  in  the  scatterplots. 
All  plots  are  valid  at  03  Z  30  December  2008. 


As  shown  in  Table  5,  DAJ  verifies  better  than  DAI  or  DA2  as  the  percent  of 
change  between  the  expected  SJ  and  actual  SJ  for  DAJ  is  larger  for  than  either  DAI  or 
DA2.  This  indicates  a  synergistic  effect  where  the  joint  assimilation  of  points  1  and  2  is 
slightly  better  than  assimilating  at  any  of  the  single  points  alone.  For  the  DAJ  run,  the 
mean  of  7  decreases  as  the  majority  of  members  move  in  a  negative  direction  (Figures  50 
and  51).  This  negative  shift  in  7  is  due  increased  shear,  decreased  stability,  and  increased 
mean  potential  temperature.  This  again  shows  that  7  continues  to  behave  in  a  consistent 
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Bulk  Richardson  manner,  and  a  decrease  in  7  corresponds  to  the  likelihood  of  turbulence 
being  increased  due  to  increased  shear  in  the  region  of  the  flow  separation.  The  increase 
in  shear  is  primarily  due  to  a  2.5  m  s'1  increase  in  the  easterlies  in  the  return  flow  under 
the  rotor.  Additionally,  the  westerlies  and  the  region  of  flow  separation  shift  slightly  to 
the  east.  This  causes  the  ensemble  mean  in  the  DAJ  run  to  increase  over  5  m  s’1  in  the 
region  just  west  of  the  Airfield  when  compared  to  the  NP  run  as  also  seen  in  the  DAI  run 
(Figure  43). 


Figure  50.  a)  Histogram  of  the  forecast  metric  (7 )  for  DAJ  96-member  WRF  ensemble  at 
03Z  30  December  2008,  b)  scatterplot  of  6  vs.  ( 92  -9X) ,  c)  scatterplot  of 

( U 2  -Ux)2  vs.  0 ,  d)  {02  -9X)  vs.  ( U2  -t/j)2  .  The  magnitude  of  the  forecast 
metric  (7)  is  shaded  in  the  scatterplots.  All  plots  are  valid  for  03Z  30 

December  2008. 
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Figure  51.  a)  Histogram  for  the  difference  in  the  forecast  metric  ( J )  between  DAJ  and 
NP  runs,  b)  scatterplot  of  the  difference  between  the  DAJ  and  NP  runs  for  6 
vs.  (<92  -<9[) ,  c)  scatterplot  of  the  difference  between  the  DAJ  and  NP  runs 

for  (U0  -t/j)2  vs.  0 ,  d)  scatterplot  of  the  difference  between  the  DAJ  and  NP 

runs  for  ( 02  -0X)  vs.  (U2  -Ux)2 .  The  magnitude  of  the  difference  between 

the  DAJ  and  NP  runs  for  the  forecast  metric  ( J)  is  shaded  in  the  scatterplots. 
All  plots  are  valid  at  03Z  30  December  2008. 


C. 


PERTURBATION  IMPACT  BEYOND  THREE  HOURS 


To  access  the  impact  each  perturbation  or  data  assimilation  had  on  the  total 
energy  of  the  inner-most  domain,  a  calculation  of  the  change  in  the  dry  total  energy  is 
performed  for  each  ensemble  run  using  the  ensemble  mean.  The  purpose  of  this 
calculation  is  to  access  the  impact  of  the  perturbation  beyond  the  three  hour  point  and  to 
gauge  the  length  of  time  ESA  may  be  useful  in  our  DSWS  scenario.  The  dry  total  energy 
is  defined  according  to  Zou  et  al.  (1997)  and  Kim  and  Jung  (2009)  as 
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where  Ed  is  the  total  dry  energy,  u  ,  v  ,  w  are  the  are  the  zonal,  meridional,  and  vertical 
wind  perturbations,  respectively;  6  is  the  potential  temperature  perturbation;  p  is  the 
pressure  perturbation;  N,  6,  p .  and  cs  arc  Brunt-Vaisala  frequency,  potential 

temperature,  density,  and  speed  of  sound,  respectively,  at  the  sea  level;  and  x,  y,  and  a 
denote  zonal,  meridional,  and  vertical  coordinate,  respectively.  Results  from  this 
calculation  are  shown  in  Figures  52,  53,  and  54. 

Figure  52  shows  the  evolution  of  total  energy  for  each  perturbed  ensemble  run  for 
the  entire  inner  domain  (domain  3)  by  forecast  hour.  Both  kinetic  and  potential  energy 
increase  in  all  ensemble  runs  for  the  first  three  hours  and  then  decrease  rapidly.  The 
increases  for  the  first  three  hours  are  noteworthy,  as  perturbations  that  are  introduced  are 
designed  utilizing  ESA  to  impact  the  forecast  at  03Z.  This  result  may  be  coincidental,  but 
is  worth  noting  for  testing  in  future  experiments  where  the  perturbation  could  be 
introduced  to  affect  forecast  metrics  at  later  forecast  hours.  These  plots  also  reveal  that 
the  full  data  assimilation  (DAI,  DA2,  and  DAJ)  introduce  perturbations  that  cause  the 
total  energy  to  be  approximately  an  order  of  magnitude  higher.  This  is  most  likely  due  to 
the  larger  Sx  introduced  by  the  assimilation  process  as  compared  to  the  perturbation 
method. 

Figure  52  also  reveals  that  the  perturbations  introduced  at  the  upstream  point 
(point  2)  impact  the  total  energy  perturbation  longer  than  the  perturbations  introduced  at 
the  downstream  (point  1).  Even  though  the  initial  total  energy  perturbation  is  smaller  in 
runs  perturbed  at  point  2,  the  decay  in  perturbation  energy  is  slower  for  runs  P2  and  DA2 
when  compared  to  PI  and  DAI.  Point  2  is  the  upstream  perturbation  and  changes 
upstream  would  take  longer  to  propagate  out  of  the  domain.  This  is  also  seen  in  Figures 
53  and  54.  Figure  53  divides  the  inner  domain  into  two  parts,  east  and  west,  and  shows 
the  total  energy  perturbation  in  these  two  regions.  The  west  calculation  is  the  western 
three-fifths  of  domain  3  or  roughly  the  mountainous  portion  of  the  domain.  The  east 
calculation  is  the  eastern  two-fifths  of  the  domain  or  roughly  the  plains  portion  of  the 
domain.  Despite  the  eastern  portion  being  one-fifth  smaller,  the  total  energy  perturbation 
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over  the  plains  is  larger  for  all  runs,  again  indicating  the  larger  impact  of  perturbations  for 
the  easterly  point  1. 

Figure  54  is  the  total  energy  perturbation  for  a  40  x  40  grid  point  box  from  the 
innermost  domain  centered  over  USAFA.  The  forecast  metric  J,  that  each  of  the 
perturbations  are  designed  to  impact  at  the  three  hour  forecast,  is  centered  in  this  box. 
These  plots  remain  consistent  with  the  previous  discussion  showing  that  the  DA  runs 
have  a  larger  energy  perturbation  and  the  perturbations  introduced  at  point  2  impact  the 
energy  perturbation  for  a  longer  period  of  time.  The  upstream  perturbations  (P2  and 
DA2)  have  a  much  longer  impact  in  this  calculation  than  any  of  the  others  presented.  In 
fact,  the  energy  perturbations  remain  larger  than  the  initial  condition  perturbations  for 
nine  hours  and  don’t  start  to  precipitously  drop  off  until  after  12  hours  into  the  forecast. 

The  total  energy  perturbation  in  the  east  for  P2  (Figure  53)  and  the  APE 
perturbation  over  USAFA  for  P2  (Figure  54)  increases  from  the  15  to  18  hour  forecast. 
This  signifies  a  potential  impact  to  the  forecast  out  to  18  hours.  This  suggests  that 
carefully  selected  and  targeted  observations  when  assimilated  into  the  analysis  can 
impact  the  ensemble  mean  forecast  for  at  least  12  hours.  The  sharp  drop  after  12  hours 
may  indicate  a  lack  of  impact  beyond  this  point.  However,  a  sharp  cold  front  with  a 
strong  synoptic  signal  begins  to  impact  the  forecast  region  at  this  time.  Although  very 
strong  and  intense,  the  DSWS  is  a  mesoscale  feature  with  relatively  small  spatial  and 
temporal  scales.  The  cold  front,  supported  by  an  upper-level  trough  and  jet  stream, 
represents  a  much  larger  spatial  and  temporal  scale  that  dominates  the  smaller  scales. 
This  cold  front  is  a  dramatic  regime  shift  and  one  expects  convergence  of  all  ensembles 
from  a  westerly  DSWS  regime  into  a  synoptic  cold  frontal  pattern  no  matter  how  big  the 
perturbations  that  are  introduced,  or  how  strong  the  DSWS. 
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Figure  52.  Kinetic,  available  potential,  and  total  energy  perturbation  for  domain  3  by 
forecast  hour  for  each  ensemble  run  from  00Z  30  December  to  00Z  3 1 
December  2008.  Units  for  the  vertical  axis  are  Joules.  The  y-axis  scales  are 

different  for  each  plot. 


89 


- East 

- West 


DAI  DA2  DAJ 


Figure  53.  Total  energy  perturbation  for  domain  3  divided  into  the  western  and  eastern 
portions  of  the  domain  by  forecast  hour  for  each  ensemble  run  from  00Z  30 
December  to  00Z  31  December  2008.  Units  for  the  vertical  axis  are  Joules. 
The  y-axis  scales  are  different  for  each  plot. 
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Figure  54.  Kinetic,  available  potential,  and  total  energy  perturbation  for  a  40  grid  point 
by  40  grid  point  box  centered  over  USAFA  from  domain  3  by  forecast  hour 
for  each  ensemble  run  from  00Z  30  December  to  00Z  31  December  2008. 
Units  for  the  vertical  axis  are  joules.  The  y-axis  scales  are  different  for  each 

plot. 


D.  METHOD  FOR  TESTING  LINEARITY  WITH  FEWER  MEMBERS 

The  method  used  in  this  research  for  testing  sensitivity  involves:  1)  creating  and 
running  a  control  ensemble  run  with  96  members  (NP  run);  2)  performing  an  ESA 
calculation  to  select  perturbation  points;  3)  perturbing  the  initial  conditions  and  then  re¬ 
running  the  entire  ensemble  (96  members)  to  determine  if  the  change  in  the  initial 
conditions  ESA  predicted  actually  changed  the  forecast  metric  in  the  same  manner.  This 
is  sometimes  referred  to  as  a  linearity  test.  This  entire  process  is  very  time  consuming  and 
expensive.  A  method  is  tested  to  try  and  simplify  this  process  and  identify  more  quickly 
if  a  perturbation  at  a  particular  location  will  impact  the  forecast  as  predicted  by  ESA  or 
not.  The  results  obtained  from  the  full  method  discussed  in  previous  sections  are  used  to 
verify  if  the  same  information  can  be  gained  by  running  just  a  few  members  rather  than 
all  96  members  in  step  3  of  this  process.  If  a  method  is  devised  to  test  the  potential 
impact  of  perturbation  points,  then  more  points  can  be  tested  in  a  timely  and  cheaper 

91 


fashion  and  only  points  that  impact  the  forecast  can  be  used  for  assimilation  (G.  Hakim, 
personal  communication,  December  14,  2012). 

The  method  tested  involves  rank  ordering  the  J  forecast  metric  by  the  member 
closest  to  the  ensemble  mean  to  the  member  farthest  from  the  ensemble  mean.  This  is 
done  for  the  NP  run  at  the  forecast  valid  time  for  J  (03Z).  This  allows  us  to  identify 
which  members  are  closest  to  the  mean  in  terms  of  the  forecast  metric  J.  Once  the  rank 
order  of  the  members  is  identified,  this  same  rank  order  of  the  members  is  applied  to  each 
of  the  perturbation  runs  (PI,  P2,  DAI,  DA2,  DAJ)  and  the  difference  between  the 
ensemble  mean  for  the  perturbation  run  and  each  individual  member  is  calculated  with 
members  remaining  in  the  original  rank  order.  This  difference  is  then  summed  and 
averaged  for  ensemble  sizes  of  1,  2,  3... etc.,  through  96  in  the  rank  order  previously 
determined.  This  calculation  reproduces  what  the  ensemble  mean  would  be  if  a  smaller 
number  of  ensembles  members  were  used  to  test  linearity.  This  calculation  is  shown 
graphically  in  Figure  55.  The  results  from  the  DAJ  run  are  omitted  as  they  are  almost 
identical  to  the  DAI  results. 

Figure  55  shows  that  a  substantial  number  of  members  are  needed  for  the 
ensemble  mean  of  smaller  sets  of  members  to  approximate  the  mean  value  of  J  for  all  96 
members  in  terms  of  the  magnitude  of  change  for  the  perturbed  runs.  The  approximate 
magnitude  of  J  for  96  members,  starts  to  occur  at  around  60  members  for  the  PI,  DAI, 
and  DAJ  runs  and  doesn’t  occur  until  greater  than  80  members  for  the  P2  and  DA2  runs. 
This  result  suggests  that  for  the  phenomena  modeled  and  the  J  chosen,  it  is  difficult  to 
capture  the  magnitude  of  the  sensitivity  of  /  in  a  perturbed  run  with  significantly  fewer 
members. 

However,  the  direction  of  change  (positive  or  negative)  of  the  entire  ensemble 

mean  is  predicted  with  far  fewer  members.  The  PI,  DAI,  and  DAJ  runs  correctly  predict 

the  direction  of  change  with  as  few  as  three  members.  The  P2  and  DA2  runs  correctly 

predict  the  right  direction  with  fewer  than  three  members,  but  then  predict  the  wrong 

direction  until  an  ensemble  of  approximately  20  members  is  obtained.  It  is  also  important 

to  note  that  P2  performed  the  worst  in  terms  of  the  linearity  tests  of  all  the  ESA  results 

(see  overview  to  this  chapter)  and  would  be  the  most  likely  perturbation  point  to  discard 
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and  replace  with  another.  These  results  suggest  that  although  smaller  number  of  ensemble 
members  may  not  be  useful  in  predicting  the  final  magnitude  of  the  change  in  J  for  the 
entire  ensemble  run,  as  few  as  three  members  may  give  insight  to  the  final  direction  of 
change  for  the  entire  set.  This  may  be  useful  in  future  work,  allowing  only  a  few  hours  to 
be  spent  testing  potential  perturbation  points  to  determine  if  the  forecast  metric  moves  in 
the  correct  direction  rather  than  days  to  weeks  running  the  entire  ensemble.  If  the  forecast 
metric  moves  in  the  wrong  direction,  results  from  this  method  suggest  that  a  new 
perturbation  point  might  be  a  wiser  choice  before  perturbing  and  running  the  full 
ensemble. 


Figure  55.  Summed  differences  of  the  forecast  metric  J  between  ensemble  runs  and  the 
NP  run  (y-axis)  for  ranked  ensemble  sizes  1-90  (x-axis). 
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V.  OBSERVATIONAL,  MODELING,  AND  FREQUENCY 
ANALYSIS  OF  FOURTEEN  SEVERE  DOWNSLOPE  WINDSTORMS 

AT  USAFA 


A.  OVERVIEW 

Despite  much  research  into  the  initiation  and  causes  of  severe  DSWS  over  many 
decades,  much  of  the  theory  is  still  based  upon  idealized  terrain  and  two-dimensional 
simulations.  While  idealized  two-dimensional  simulations  give  great  insight  into  what 
atmospheric  parameters  play  a  role  in  the  development  of  these  storms,  they  lack  full 
explanations  of  how  these  theoretical  principles  play  out  in  more  realistic  complex  three- 
dimensional  flows.  The  goal  of  the  work  done  in  this  chapter  is  to  determine  whether  the 
conditions  and  atmospheric  features  that  caused  the  severe  DSWS  on  30  December  2008, 
and  evaluated  using  ESA,  are  unique  or  not  when  put  in  context  of  other  DSWS.  The 
research  aims  to  model,  analyze,  and  classify  the  strongest  DSWS  recorded  by  the 
HWAS  sensors  at  USAFA.  This  helps  determine  which  of  the  downslope  initiation 
factors  hypothesized  in  the  literature  are  most  important  at  this  location,  which 
classification  of  storms  are  the  strongest,  and  which  have  winds  that  propagate  the 
farthest  from  the  lee  of  the  mountain  slope. 

B.  BACKGROUND 

An  in-depth  review  of  downslope  winds,  and  overview  of  the  historical  research, 
is  given  in  Chapter  II.  Durran  (2003)  summarizes  and  states  that,  “Significant  downslope 
winds  have  been  observed  to  develop  in  three  basic  situations:  (1)  when  a  standing 
mountain  wave  in  a  deep  cross-mountain  flow  achieves  sufficient  amplitude  to  overturn 
and  breakdown  at  some  level  in  the  troposphere,  (2)  when  standing  mountain  waves 
break  and  dissipate  at  a  critical  level  in  a  shallow  cross-mountain  flow,  and  (3)  when 
there  is  sufficient  static  stability  near  mountain-top  level  in  the  cross-mountain  flow  to 
create  high  downslope  winds  even  without  wave  breaking.”  For  the  purposes  of  this 
research,  type  1  from  above  will  be  referred  to  as  wave  breaking  (WB).  Type  2  events  are 
most  common  during  a  bora  and  are  not  examined,  as  they  are  not  found  to  be 

responsible  for  the  strongest  DSWS  at  USAFA.  Type  3  storms,  or  the  non-breaking 
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cases,  often  occur  when  there  is  a  strong  stable  layer  near  mountaintop  level,  below  a 
much  deeper  layer  of  relatively  weak  upper-tropo spheric  stability  (Durran  1986).  These 
non-breaking,  stability  layered,  less  vertical  type  3  cases  are  referred  to  as  non-wave 
breaking  (NWB)  for  this  research. 

Prototypical  examples  of  vertically  propagating  breaking  waves  and  stability 
layered  non-breaking  trapped  waves  are  given  in  Durran  2003  for  idealized  isolated 
terrain  and  a  two-dimensional  model.  The  strong  downslope  wind  caused  by  a  vertically 
propagating  breaking  wave  is  characterized  by  a  positive  horizontal  velocity  perturbation 
of  greater  than  18  ms'1  (Figure  5  Durran  2003)  and  a  maximum  surface  wind  speed  of 
greater  than  28  m  s"1  is  noted.  The  prototypical  example  of  a  non-breaking  trapped  wave 
(Figure  6  Durran  2003)  shows  non-breaking  waves  lead  to  strong  winds,  but  significantly 
weaker  winds  than  the  wave  breaking  case,  with  the  positive  horizontal  velocity 
perturbation  of  only  12  ms'1  for  this  simulation  along  the  lee  slope. 

These  idealized  cases  have  led  to  the  popular  notion  that  the  strongest  DSWS  are 
almost  always  associated  with  vertically  propagating  wave  breaking  events.  Evidence  can 
be  found  in  educational  tools  such  as  the  popular  COMET  modules  produced  by 
University  Corporation  for  Atmospheric  Research  (UCAR).  The  COMET  module 
entitled  “Mountain  Waves  and  Downslope  Winds,”  states  that  strong  downslope  wind 
cases  are  usually  associated  with  waves  breaking  aloft  and  that  the  strongest  downslope 
winds  develop  directly  beneath  the  region  of  flow  reversal  at  the  critical  level  in  mature 
vertically  propagating  waves  (METED/COMET  cited  2013).  While  it  is  true  that 
vertically  propagating  waves  with  a  critical  level  and  wave  breaking  often  produce  strong 
winds,  the  research  presented  in  this  chapter  suggests  that  in  the  highly  irregular  and 
three-dimensional  terrain  upstream  of  USAFA,  the  strongest  events  recorded  by  the 
HWAS  sensors  are  NWB.  Winds  from  these  NWB  events  on  average  extend  farther  to 
the  east  and  away  from  the  lee  slope  than  vertically  propagating  WB  events. 
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C.  DATA  AND  METHODOLOGY 

The  severe  DSWS  that  occurred  on  30  December  2008  at  USAFA,  that  is 
modeled  and  studied  by  the  WRF  ensemble  runs  for  the  ESA  research  presented  in 
previous  chapters,  is  the  strongest  and  has  the  largest  extent  of  strong  winds  ever 
recorded  by  the  HWAS  network.  As  noted  in  Chapter  III,  this  storm  is  caused  by  a  NWB 
trapped  wave  upstream  of  the  Front  Range  over  the  higher  terrain  of  the  central  Colorado 
Rockies.  To  understand,  whether  or  not  this  type  of  NWB  event  is  commonplace  at 
USAFA,  thirteen  additional  severe  DSWS  are  identified  using  the  HWAS  database  and 
deterministically  modeled  using  WRF.  The  fourteen  storms  in  total  are  chosen  based 
upon  the  strongest  winds  observed  at  the  various  HWAS  sensors.  In  some  cases  severe 
winds  are  observed,  but  only  at  a  handful  of  sensors. 


Date 

#  of 

HWAS 

Severe 

Max  Gust 

Rampart 

(kt) 

Max  Gust 

Community 
Center  (kt) 

Max  Gust 

Airfield  (kt) 

Max  Gust 

Pine 

Creek  (kt) 

Strength 

Score 

Type 

Upstream 

Wavelength 

(km) 

Max  Wind 

Model  (kt) 

3/6/2004 

6 

73 

58 

45 

41 

217 

NWB 

48.78 

82 

11/3/2005 

12 

71 

70 

63 

63 

267 

NWB 

74 

4/15/2006 

9 

73 

69 

62 

46 

250 

WB 

74 

6/7/2007 

11 

60 

63 

58 

51 

232 

WB 

15.50 

66 

10/2/2007 

4 

68 

63 

34 

37 

202 

47 

1/5/2008 

11 

74 

65 

54 

57 

250 

74 

12/18-19/2008 

5 

72 

57 

47 

41 

WB 

26.25 

74 

12/30/2008 

12 

79 

68 

64 

57 

268 

NWB 

26.18 

62 

12/20/2010 

10 

62 

56 

57 

53 

228 

NWB 

42.74 

66 

11/12/2011 

10 

68 

68 

60 

52 

248 

NWB 

26.25 

66 

11/13/2011 

4 

77 

52 

49 

43 

221 

WB 

16.26 

86 

1/19/2012 

73 

66 

53 

54 

246 

NWB 

48.78 

82 

2/21/2012 

7 

61 

60 

60 

37 

218 

WB 

14.84 

54 

2/23/2012 

11 

66 

60 

60 

51 

237 

WB 

42.74 

78 

Table  7.  Table  summarizing  fourteen  severe  DSWS  analyzed  and  modeled  at 
USAFA.  Six  non-wave  breaking  (NWB)  and  eight  wave  breaking  (WB) 

storms  are  shown. 


Table  7  includes  a  list  of  these  fourteen  severe  DSWS,  and  parameters  describing 
these  storms.  In  total,  only  two  storms  have  ever  produced  severe  winds  at  all  twleve 
HWAS  sensors  and  only  eight  storms  have  been  observed  at  10  or  more  sensors.  Each 


97 


0  30  60  90  120  150  180 

A  -  Smmtch  Mountains  D  -  South  P.rt  X-TflLsl  COIUpOllCllT  | 111  «  1 


B  -  Arkansas  River  VbSey 
C -  Mosquito  Range 


C  -  Divide 

F-  Front/ Ram  part  Range 


■JiiiiMiiirmi 


-4  0  1  8  12  16  20  21  28  32  36 


1 

•-  Rampart  O  " 

•-  Community  Center  #-  Pine  Creek 


x-wmd  component  (m  s-1) 

Mill!  M"il  IMI 

-4  0  4  8  12  16  20  24  28  32  36 


Figure  56.  Cross  section  of  .r-wind  component  (shaded)  and  potential  temperature 
(contoured)  from  west  to  east  across  Colorado  (top)  at  the  gridline  of  the 
Rampart  sensor  from  the  WRF  deterministic  run  valid  23Z  18  December 
2008.  The  black  box  indicates  the  area  shown  in  the  cross  section  in  the 
bottom  part  of  the  figure.  The  horizontal  axis  is  scaled  to  the  number  of  grid 
points  from  the  innermost  domain  with  a  resolution  of  1.33  km.  Approximate 
locations  of  four  HWAS  sensors  are  marked  by  circles. 


98 


of  these  eight  storms  are  included  in  the  set,  as  well  as  six  other  storms  notable  for  their 
extreme  winds  at  one  or  more  of  the  HWAS  sensors. 

Single  deterministic  WRF-ARW  runs  are  performed  for  each  of  these  severe 
DSWS.  The  truth  run  from  the  30  December  2008  storm  is  selected  as  the  model  of 
choice  for  that  event.  The  configuration  of  the  truth  run  (domains,  resolution,  physics 
packages,  etc.)  is  used  for  the  deterministic  WRF  runs  (c.f.  Chapter  III,  Figure  20,  Table 
1).  The  forecast  is  output  at  one  hour  intervals.  Cross  sections  of  x-wind  component  and 
potential  temperature  are  plotted  for  each  event  at  every  forecast  hour  along  the  solid 
black  line  given  in  Figure  16.  In  addition,  a  subdomain  is  plotted  at  each  forecast  hour, 
showing  more  detail  in  the  area  given  by  the  black  box  in  Figure  56.  Each  event  is 
analyzed  hour  by  hour  to  determine  the  forecast  time  when  the  model  predicts  the 
strongest  x-wind  component  in  this  box.  The  time  and  speed  of  the  forecast  maximum 
winds  around  USAFA  is  recorded.  The  speed  of  the  maximum  wind  for  each  event  as 
recorded  by  HWAS  at  the  Rampart,  Community  Center,  Airfield,  and  Pine  Creek  sensors, 
and  forecast  by  the  model  in  the  zoomed  in  cross  sections  is  given  in  Table  7. 

These  cross  sections  are  analyzed  and  each  storm  is  classified  into  one  of  two 
general  categories:  non-wave  breaking  (NWB)  or  wave  breaking  (WB)  events.  The 
primary  criteria  used  to  determine  storm  category  is  the  existence  of  middle  to  upper- 
tropospheric  self-induced  critical  layers,  as  evidenced  by  easterlies  or  negative  x-wind 
components  in  the  cross  sections.  If  a  self-induced  critical  layer  is  evident,  suggesting  a 
vertically  propagating  wave  with  wave-breaking,  then  the  storm  is  categorized  as  wave 
breaking  (WB).  If  no  self-induced  critical  layer  is  found,  then  the  storm  is  categorized  as 
non-wave  breaking  (NWB).  The  categorization  of  each  storm  is  given  in  Table  7.  Six  of 
the  storms  were  categorized  as  non-wave  breaking  and  eight  of  the  storms  were 
categorized  as  wave-breaking.  As  discussed,  the  most  notable  non-wave  breaking  storm 
was  the  30  December  2008  event,  used  for  the  ESA  portion  of  this  research  and  depicted 
by  the  truth  run  shown  in  Appendix  B .  An  example  of  a  wave  breaking  storm  is  shown  in 
Figure  56  and  cross  sections  at  the  time  of  maximum  winds  are  given  in  Appendix  C. 

The  criterion  used  to  identify  NWB  events  is  the  lack  of  evidence  in  the  model  of 

a  self-induced  critical  level.  While  the  non-existence  of  a  vertically  propagating  breaking 
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wave  does  not  automatically  determine  that  a  trapped  wave  exists,  evidence  suggests  that 
the  majority  of  the  non-wave  breaking  events  studied  for  this  research,  exhibit  trapped 
behavior  or  properties.  Note  that  the  term  “trapped”  used  here  refers  to  trapping  of  wave 
energy  upstream  of  the  lee  slopes  and  location  of  the  severe  DSWS,  not  trapped  waves  in 
the  lee  of  the  terrain.  Trapped  lee  waves,  rotor,  flow  separation,  and  flow  reversal  at  the 
surface  under  the  rotor  are  observed  and  suggested  by  both  the  surface  HWAS  data  and 
the  model  simulations  downstream  of  the  Front  Range  for  a  period  of  time  in  almost 
every  event,  but  are  not  used  as  a  criteria  for  determining  the  type  of  storm  causing  the 
downslope  winds. 

D.  RESULTS 

1.  Strength  Score  Analysis 

A  strength  score  is  developed  for  downslope  winds  observed  at  USAFA.  This 
score  is  the  sum  of  the  peak  wind  gusts  observed  during  each  storm  by  the  Rampart, 
Community  Center,  Airfield,  and  Pine  Creek  HWAS  sensors.  These  sensors  are  chosen 
as  they  are  all  roughly  at  the  same  latitude  and  extend  in  a  line  from  west  to  east  allowing 
for  the  assessment  of  the  strength  and  distance  from  the  lee  slope  the  strong  winds 
propagate.  The  exact  location  of  these  sensors  is  given  in  Figure  63  and  Table  8  in 
Appendix  A.  The  maximum  gust  values  for  each  storm,  at  each  sensor,  and  the  summed 
strength  score  is  given  in  Table  7. 

The  strength  score  reveals  that  the  strongest  winds  recorded  at  each  of  these  four 
sensors  all  occurred  during  a  trapped  event.  The  highest  recorded  gusts  at  the  Rampart 
and  Airfield  sensors  (79  and  64  kt)  occurred  during  the  30  December  2008  NWB  event 
and  the  highest  recorded  gusts  at  the  Community  Center  and  Pine  Creek  sensors  (70  and 
63  kt)  occurred  during  the  03  November  2005  NWB  event.  These  two  events  rank  a  close 
first  and  second  in  overall  strength  score  with  values  of  268  and  267,  respectively.  Note 
that  these  two  NWB  events  are  the  only  two  DSWS  to  record  severe  wind  gusts  (greater 
than  50  kt)  at  all  twelve  HWAS  sensors.  Overall,  four  of  the  highest  six  strength  scores 
are  produced  by  NWB  events  rather  than  vertically  propagating  WB  events.  In  addition, 
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the  average  strength  score  for  NWB  events  is  much  higher  than  vertically  propagating 
WB  events  with  a  score  of  246  compared  to  228. 


Figure  57.  Composite  vertical  soundings  by  model  level  at  Grand  Junction,  CO  for  all 
non-breaking  events  (top)  and  breaking  events  (bottom).  Valid  three  hours 
prior  to  strongest  winds  at  USAFA.  Shown  is  Potential  temperature  (K)  in 
plots  a  and  e,  x-wind  component  (m  s’1)  in  plots  b  and  f,  Brunt-Vaisala 
frequency  in  plots  c  and  g,  and  Scorer  parameter  in  plots  d  and  h. 

2.  Composite  Upstream  Soundings 

Composite  soundings  of  all  NWB  and  WB  events  are  calculated  at  the  upstream 
sounding  location  of  Grand  Junction  and  are  shown  in  Figure  57.  Composite  soundings 
are  valid  three  hours  prior  to  the  time  of  the  strongest  modeled  winds  at  USAFA.  These 
soundings  show  that  NWB  events  are  characterized  by  stability  layering  as  revealed  by 
potential  temperature  and  the  Brunt-Vaisala  frequency.  In  the  NWB  composites,  a  stable 
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layer  is  present  as  revealed  by  a  relative  maximum  of  Brunt- Vaisala  frequency  at  model 
levels  12  to  18,  which  corresponds  to  an  altitude  of  about  3500  -  5500  m  MSL.  In 
addition,  the  NWB  events  have  stronger  forward  shear  that  inhibits  wave  breaking  due  to 
a  much  sharper  increase  in  westerly  winds  with  height,  with  an  average  maximum  wind 
of  46  m  s"1  centered  at  model  level  35  (11  km  MSL).  These  factors  lead  to  a  strong 
increase  with  height  of  the  Scorer  parameter  in  the  lower  and  middle  troposphere  (Figure 
57d)  that  is  commonly  associated  with  trapped  waves.  In  contrast,  the  WB  events  have 
only  one  pronounced  stable  layer  near  the  ground,  winds  that  increase  much  slower  with 
height,  and  a  much  weaker  wind  maximum  of  only  28  m  s'1  at  model  level  34  (10.7  km 
MSL). 


Vertical  Velocity  -  Level  10  from  (200.178)  to  (350.178) 


Figure  58.  Upstream  z-wind  component  at  the  gridline  of  the  Rampart  sensor  during  time 
of  maximum  winds  at  USAFA  from  model  level  10  for  the  truth  run  valid 

12Z  December  2008. 
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AmpMude  Spectrum  of  Vertical  Velocity 
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Figure  59.  Frequency  spectrum  of  upstream  z-wind  component  during  time  of  maximum 
winds  at  USAFA  from  model  level  10  for  the  truth  run  valid  12Z  December 

2008. 


3.  Fourier  Analysis 

To  access  the  upstream  wavelength  properties  of  NWB  and  WB  events,  and  to 
compare  these  results  to  the  dominant  upstream  terrain  frequencies,  a  Fourier  analysis  is 
conducted  of  the  upstream  waves  and  terrain.  The  purpose  of  this  analysis  is  to  determine 
what  role  the  wavelength  of  the  terrain  has  in  determining  the  type  of  gravity  wave 
causing  the  DSWS  and  if  wavelength  varies  between  NWB  and  WB  events. 

To  determine  the  dominant  frequencies  and  wavelength  of  the  upstream  waves,  a 
Fourier  analysis  of  the  upstream  wave  structures  is  performed  on  all  fourteen  wind 
events,  by  conducting  a  Fast  Fourier  Transform  (FFT)  of  the  z-wind  component  at  model 
level  10  (approximately  1400  m  AGL)  from  grid  points  200  to  350  of  domain  3  at  the 
gridline  of  Rampart  sensor  at  the  hour  of  maximum  winds  in  the  model  at  USAFA.  Since 
the  data  varies  as  a  function  of  distance  (km)  rather  than  time,  the  FFT  transforms  the 
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data  into  wave  number  space,  giving  frequencies  with  units  of  km"1  or  the  number  of 
wave  cycles  per  km.  The  wavelength  with  the  largest  amplitude  in  the  frequency  domain 
is  given  in  Table  7  for  each  event.  The  "-wind  component  at  the  time  of  maximum  winds 
at  USAFA  for  the  trapped  event  on  30  December  2008  from  the  truth  run  is  given  in 
Figure  58  and  the  upstream  frequency  domain  for  this  event  is  given  in  Figure  59. 

A  Fourier  analysis  of  the  upstream  terrain  height,  as  resolved  by  the  model,  is  also 
conducted  by  performing  a  FFT  at  the  gridline  of  the  Rampart  sensor  and  for  the  same 
grid  points  as  the  wind  analysis.  A  method  is  used  that  follows  Young  and  Pielke  (1983), 
where  a  linear  trend  is  calculated  in  the  terrain  by  the  method  of  least  squares  and  then 
subtracted  from  the  terrain  height.  The  resulting  series  is  transformed  into  wavenumber 
space  using  a  FFT  and  frequencies  are  given  as  the  number  of  cycles  per  km.  Figure  60 
shows  a  plot  of  the  frequency  domain  of  the  upstream  terrain.  The  wavelength  of  each  of 
the  fourteen  storms  analyzed  is  depicted  on  the  frequency  plot  as  either  a  for  a  NWB 
event  or  as  a  “#”  for  a  WB  event.  For  the  strongest  severe  DSWS  there  are  three  distinct 
modes  or  frequencies  that  account  for  all  of  these  storms,  and  each  of  the  three  modes 
correspond  to  dominant  frequencies  in  the  terrain  spectrum.  The  three  active  modes 
correspond  to  a  wavelength  of  42.5  km,  26.2  km,  and  14.8  km.  The  trapped  events  are 
confined  to  the  longer  two  modes,  and  the  vertical  events  are  mostly  confined  to  the 
shorter  two  modes. 

These  results  are  in  contrast  to  previous  work,  showing  that  the  wavelength  of 
breaking  events  at  Boulder,  Colorado  have  a  longer  wavelength  than  the  trapped  non¬ 
breaking  waves  downstream  (Klemp  and  Lilly  1975).  Analysis  here  shows  that  the 
wavelength  of  the  NWB  events  causing  DSWS  is  on  average  longer  than  the  WB  events. 
Figure  61  reinforces  this  conclusion  with  the  terrain  frequency  domain  plotted  against  the 
composite  NWB  and  WB  frequency  domains.  This  plot  also  shows  the  alignment  of  the 
composite  NWB  storms  with  the  two  longer  wavelength  active  terrain  modes,  and  the 
alignment  of  the  composite  WB  storms  with  the  two  shorter  wavelength  active  terrain 
modes.  This  plot  also  reveals  that  there  are  dominant  terrain  modes  such  as  the  one  at  114 
km  that  do  not  correspond  to  active  modes  for  severe  DSWS  in  the  atmosphere. 
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Figure  60.  Frequency  spectrum  of  upstream  terrain  overlaid  with  wavelength  of  the 

fourteen  severe  DSWS  as  given  in  Table  7.  A  denotes  the  wavelength  of 
a  non-wave  breaking  event  and  a  “#”  denotes  the  wavelength  of  a  wave 

breaking  event. 
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Figure  61.  Frequency  spectrum  of  upstream  terrain  (green)  and  composite  frequency 
spectrums  for  the  non-wave  breaking  (blue)  and  wave  breaking  (red)  severe 

DSWS. 


E.  SUMMARY  AND  DISCUSSION 

This  chapter  shows  the  strongest  DSWS,  both  in  magnitude  and  eastward 
propagation  away  from  the  mountains,  are  caused  by  NWB  events  far  more  often  than 
WB  events  at  USAFA.  The  majority  of  the  NWB  events  exhibit  trapped  behavior  or 
properties.  This  suggests  that  trapped  waves,  perturbed  by  the  upstream  terrain  over 
central  Colorado,  propagate  downstream  across  the  elevated  terrain  west  of  the  Rampart 
Range  and  play  a  significant  role  in  the  formation  of  severe  DSWS  at  USAFA.  This 
highlights  the  theory  adopted  by  many  that  the  strongest  DSWS  are  caused  by  vertically 
propagating  breaking  waves.  Much  of  this  theory  is  derived  from  two-dimensional 
regular  terrain  and  idealized  simulations.  Our  research  suggests  that  strong  downslope 
winds  in  complex  terrain  may  be  caused  by  upstream  wave  trapping  more  often  than 
previously  thought,  and  that  broader  climatological  classification  studies  at  more 
locations  may  help  clarify  severe  DSWS  dynamics. 
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The  30  December  2008  severe  DSWS  that  was  caused  by  a  NWB  trapped  wave 
was  not  only  isolated  to  the  irregular  terrain  around  USAFA.  As  noted  in  Chapter  II,  this 
windstorm  was  widespread  causing  more  than  $7  million  in  damage  in  Colorado,  making 
it  the  fourth  most  expensive  windstorm  in  Colorado’s  history.  Cross  sections  taken  from 
the  truth  run  at  Boulder,  CO  reveal  that  a  NWB  trapped  wave  is  also  responsible  for  the 
high  winds  and  damage  at  that  location  which  is  characterized  by  more  regular  terrain 
and  a  continuous  decent  form  the  Continental  Divide  to  the  Eastern  Plains.  The  highest 
wind  gusts  recorded  in  Boulder  from  this  upstream  trapped  wave  were  96  kt  (Earth 
System  Research  Laboratory,  cited  2013).  This  suggests  that  extremely  strong  NWB 
trapped  events  may  play  a  more  significant  role  at  locations  of  both  irregular  and  more 
idealized  terrain. 

This  study  highlights  the  differences  in  upstream  atmospheric  soundings  that 
cause  these  two  different  types  of  storms  as  shown  in  the  composite  soundings  from 
Grand  Junction.  The  results  from  this  study  and  the  ability  to  determine  a  priori  using 
observed  and  forecast  upstream  soundings  should  aid  forecasters  at  USAFA  and  similar 
locations  to  develop  forecast  rules  of  thumbs  that  distinguish  between  NWB  and  WB 
events.  This  information  can  be  used  to  help  determine  the  likely  hood  that  locations 
away  from  the  base  of  the  mountains  will  experience  extreme  threshold  events. 

The  results  from  the  Fourier  analysis  of  these  severe  DSWS  show  that  only 
certain  wavelengths  are  excited  in  the  strongest  downslope  wind  events  at  USAFA.  It 
remains  to  be  seen  if  only  strong  events  occur  at  these  wavelengths  or  if  weaker  events 
also  occur  at  these  excitable  modes.  The  question  also  remains  about  whether  the  terrain 
forces  the  atmospheric  flow  into  these  preferred  wavelengths  or  if  the  forcing  comes  from 
the  atmosphere,  such  that  when  certain  parameters  like  stability,  wind  shear,  and 
tropopause  height  come  together,  the  atmospheric  dynamics  alone  determine  these 
preferred  modes.  It  is  possible  that  if  the  dynamics  of  the  atmosphere  produce  a  wave  at  a 
mode  that  is  not  preferred,  when  the  flow  interacts  with  the  terrain,  the  flow  is  sheared,  or 
there  is  destructive  interference  preventing  the  most  extreme  downslope  wind  events 
from  occurring.  A  more  in-depth  analysis  of  both  strong  and  weak  storms  could  lead  to 
the  development  of  forecast  tools  based  upon  verification  of  NWP  forecasts  that  use  the 
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wavelength  or  frequency  as  predicted  by  the  model  to  determine  the  probability  that  a 
strong  threshold  level  event  that  propagates  away  from  the  mountain  slope  is  more  likely 
to  occur. 

The  goal  of  the  analysis  presented  in  this  chapter  is  to  determine  whether  the 
conditions  and  atmospheric  features  that  caused  the  severe  DSWS  on  30  December  2008, 
and  evaluated  using  ESA,  are  unique  or  not  when  put  in  context  of  other  DSWS.  The 
results  presented  in  this  chapter  link  the  atmospheric  structures  found  in  NWB  events  to 
the  structures  successfully  identified  by  ESA.  ESA  identified  the  sensitivity  of  a  NWB 
DSWS  at  USAFA  to  a  strong  stable  layer  at  model  level  13  (2200  m  AGL)  at  both 
upstream  and  downstream  locations.  The  vertical  profile  of  the  ensemble  mean  upstream 
sensitivity  location  (Figure  26  bottom)  is  very  close  in  proximity  and  structure  to  the 
composite  sounding  at  Grand  Junction,  CO,  produced  from  six  NWB  DSWS  (Figure  57). 
The  composite  sounding  (Figure  62)  at  the  downstream  sensitivity  location  (point  1)  from 
the  six  NWB  events  shows  a  strong  stable  layer  at  model  level  15  (2800  m  AGE)  as 
revealed  by  potential  temperature  and  Brunt- Vaisala  frequency.  Results  show  that  ESA 
successfully  identified  an  atmospheric  structure  that  is  common  in  the  majority  of  the 
strongest  NWB  DSWS  at  USAFA.  This  suggests  that  targeting  observations  at  these 
locations  prior  to  the  onset  of  a  non-breaking  severe  DSWS  would  impact  and  potentially 
improve  forecasts  of  these  phenomena  in  numerical  models. 
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Figure  62.  Composite  vertical  soundings  by  model  level  at  downstream  sensitivity 

location  (point  1)  for  six  non-breaking  severe  DSWS.  Valid  three  hours  prior 
to  strongest  winds  at  USAFA.  Shown  is  a)  potential  temperature  (K),  b)  x- 
wind  component  (m  s'1),  c)  Brunt-Vaisala  frequency,  and  d)  Scorer 

parameter. 
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VI.  CONCLUSIONS  AND  FUTURE  WORK 


A.  SUMMARY 

We  present  as  primary  and  secondary  goals:  1)  to  determine  if  a  mesoscale 
ensemble  and  ESA  is  a  useful  analysis  tool  to  determine  the  sensitivity  of  severe  DSWS 
in  complex  terrain  to  initial  conditions  and  assimilated  observations,  and  2)  to  provide 
context  for  the  ESA  case  study  with  analysis  of  thirteen  additional  severe  DSWS.  ESA 
was  utilized  to  quantify  the  sensitivity  of  a  severe  DSWS  with  respect  to  leeside 
conditioning  in  contrast  to  upstream  synoptic  conditions.  This  study  aimed  to  determine 
if  sensitivities  can  be  used  to  reduce  forecast  uncertainty  for  a  forecast  metric  utilizing 
ESA  and  a  targeted  network  design.  This  study  is  motivated  by  the  larger  question  of 
whether  ESA  can  be  useful  tool  in  fine-scale,  highly  non-linear,  high  impact  events. 

These  goals  were  accomplished  by  identifying  and  studying  a  severe  DSWS  on  30 
December  2008  at  USAFA,  lee  of  the  Rampart  Range,  just  north  of  Colorado  Springs, 
Colorado.  A  mesonet  of  twelve  surface  observations  known  as  the  High  Wind  Alert 
System,  or  HWAS,  was  utilized  in  identifying  and  understanding  the  spatial  and  temporal 
impacts  of  this  strong  severe  DSWS  that  buffeted  the  local  area  with  wind  gusts  topping 
out  at  79  kt.  A  series  of  over  100  WRF-ARW  single  deterministic  runs  were  performed  to 
test  the  time  step,  vertical  resolution,  and  physics’  packages  best  suited  to  resolve  this  and 
other  DSWS  at  this  location.  A  final  WRF  configuration  was  set  and  used  in  conjunction 
with  an  Observing  System  Simulation  Experiment  (OSSE)  to  create  a  five  day  truth  run 
for  this  event.  Boundary  conditions  and  soil  conditions  were  updated  every  three  hours 
with  CFSR  data  for  this  truth  run. 

Next,  the  WRF  ensemble  system  was  created  and  conditioned  by  assimilating 
synthetic  observations  from  the  truth  run  that  mimicked  a  network  of  rawinsonde  and 
ASOS  locations  every  three  hours.  An  ensemble  referred  to  as  the  NP  run  was  created 
with  96  members  and  run  from  00Z  30  December  to  00Z  3 1  December  for  the  duration  of 
this  event.  Ensemble  Sensitivity  Analysis  was  performed  with  various  forecast  metrics  on 
the  NP  run.  A  forecast  metric  (7)  that  is  a  modified  Bulk  Richardson  number  was  settled 
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upon  and  two  points  of  strong  sensitivity  in  the  initial  conditions  (one  upstream  and  one 
downstream)  were  identified  from  ESA.  Five  separate  ensemble  runs  were  generated  (PI, 
P2,  DAI,  DA2,  DAJ).  The  first  two  (PI  and  P2)  introduced  separately  a  perturbation  of 
potential  temperature  from  the  truth  run  at  the  downstream  and  upstream  locations  and 
the  perturbation  was  spread  utilizing  linear  regression  of  the  analysis  as  sampled  by  the 
ensemble.  The  final  three  ensemble  runs  (DAI,  DA2,  DAJ)  were  generated  by 
assimilating  a  temperature  observation  from  the  truth  run  at  the  downstream  location, 
upstream  location,  and  both  locations  simultaneously  and  performing  full  data 
assimilation  using  an  Ensemble  Kalman  filter.  The  results  from  these  experiments  were 
used  to  access  our  research  objectives. 

B.  OUTCOMES  AND  BROADER  IMPLICATIONS 

The  most  important  outcomes  of  this  research  is  demonstration  of  the  utility  of 
ESA  and  its  capability  of  quantifying  severe  DSWS  sensitivity.  In  varying  degrees,  ESA 
showed  skill  in  identifying  regions  of  initial  condition  sensitivity  to  a  forecast  metric  J. 
Every  experiment  moved  the  ensemble  mean  for  J  in  the  correct  direction  as  predicted  by 
the  sensitivity  ( clJ/clx )  given  in  Table  2  and  shown  in  Figure  24.  In  general,  experiments 
at  point  1  (downstream)  had  a  larger  actual  5/  and  impact  on  the  forecast  than 
experiments  at  point  2  (upstream)  for  both  the  perturbation  experiment  and  data 
assimilation  experiment.  This  was  partially  due  to  the  larger  5x  introduced  for  DAI,  and 
also  showed  that  the  perturbation  associated  with  the  downstream  condition  is  more 
linearly  related  to  the  downslope  winds  three  hours  later. 

The  PI  ensemble  run  outperformed  the  other  model  runs  by  moving  J  over  48% 
toward  the  expected  value  as  predicted  by  the  ESA  performed  on  the  NP  run  as  shown  in 
Table  5.  While  the  perturbations  introduced  were  very  small,  even  these  small 
perturbations,  identified  by  ESA,  impacted  the  forecast  metric  J  reasonably  well  showing 
the  strong  dependence  of  mountain  waves  on  the  initial  conditions.  The  better 
performance  by  DAJ  indicated  a  synergistic  effect  where  the  joint  assimilation  of  points  1 
and  2  is  slightly  better  than  assimilating  at  any  of  the  single  points  alone. 
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Another  outcome  of  this  research  is  demonstrating  that  ESA  can  identify  physical 
and  dynamical  structures  in  the  atmosphere  that  play  a  role  in  the  formation  of  a  DSWS. 
ESA  was  able  to  identify  traditional  upstream  sounding  parameters,  such  as  the  strength 
of  the  near  mountain  top  stable  layer,  as  playing  a  role  in  the  strength  and  structure  of  a 
DSWS.  Most  importantly,  ESA  shed  new  light  on  the  impact  of  downstream  structures, 
or  leeside  conditioning,  that  play  a  role  in  the  formation  of  flow  separation,  leeside  rotor, 
and  convergence  between  the  rotor’s  return  flow  and  the  downslope  winds.  In  fact, 
observation  assimilation  that  adjusted  the  strength  and  characteristics  of  the  downstream 
air  mass  and  stable  layers  had  more  impact  on  the  forecast  than  assimilating  upstream 
observations.  Intuitively,  this  relationship  makes  sense.  While  a  DSWS  can  only  occur  if 
the  upstream  terrain,  flow  characteristics,  and  forcing  is  adequate,  once  a  wave  is 
generated,  the  characteristics  of  that  wave  and  flow  will  be  modified  by  the  environment 
it  is  flowing  into.  The  results  suggest  that  leeside  conditioning  is  important  if  trying  to 
determine  the  extent  of  westerly  winds  down  the  slope,  the  point  of  flow  separation,  and 
the  strength,  location,  and  structure  of  leeside  rotors. 

Another  outcome  of  this  research  is  highlighting  the  usefulness  of  a  modified 
Bulk  Richardson  number  in  a  less  than  traditional  application.  Traditionally,  Bulk 
Richardson  numbers  are  used  on  scales  smaller  than  studied  in  this  research  to  classify 
turbulent  vertical  motions  and  buoyancy  forces.  While  there  is  definitely  turbulence  and 
vertical  motion  in  the  dynamic  DSWS  and  lee  wave  rotor  scenario,  the  modified  Bulk 
Richardson  number  developed  for  this  study  proved  well  suited  to  capture  the  relevant 
physics  and  structures  in  the  convergence  zone  between  the  downslope  winds  and  leeside 
rotor  in  the  vicinity  of  the  flow  separation.  The  success  of  ESA  in  this  research  can  be 
directly  attributed  to  the  Bulk  Richardson  number  being  an  appropriate  metric  to  sample 
these  features. 

The  extent  to  which  sensitivities  can  be  used  to  reduce  forecast  uncertainty  is  still 
unclear.  This  study  showed  that  sensitivities  from  ESA  certainly  can  affect  the  spread  of 
an  ensemble  as  demonstrated  by  changes  in  the  standard  deviation  of  the  distribution. 
While  this  shows  sensitivity  to  the  initial  conditions,  it  does  not  necessarily  show  a 
reduction  in  forecast  uncertainty.  The  ensemble  runs  that  verified  best  in  terms  of  moving 

113 


J  toward  what  ESA  predicted  had  an  increase  in  the  standard  deviation  indicating 
possible  increased  uncertainty  in  the  forecast.  In  general,  perturbations  introduced  and 
observations  assimilated  upstream  decreased  the  spread  of  the  forecast  metric  increasing 
confidence  in  the  forecast.  Perturbations  introduced  or  observations  assimilated 
downstream  increased  the  spread  of  the  forecast  metric  decreasing  confidence  in  the 
forecast.  An  explanation  is  proposed  and  supported  by  evidence  that  shows  advection  of 
regions  of  decreased  spread  downstream  from  the  assimilation  points  leading  to 
decreased  spread  from  the  upstream  point  and  increased  spread  from  the  downstream 
point. 

The  modeling,  composite,  and  Fourier  analysis  performed  on  fourteen  severe 
DSWS  showed  that  upstream  non-breaking  waves,  that  exhibit  trapped  behavior  and 
properties,  play  a  larger  role  in  severe  DSWS  at  USAFA  than  previously  thought.  This 
study  showed  that  non-breaking  waves  are  responsible  for  the  strongest  winds  recorded 
by  HWAS  and  a  majority  of  the  strongest  DSWS  that  propagate  away  from  the  mountain 
slope  and  out  over  the  plains  despite  the  fact  that  severe  wave  breaking  DSWS  occur 
more  frequently.  The  Fourier  analysis  revealed  that  all  fourteen  severe  storms  fall  into 
one  of  three  frequency  modes  that  correspond  to  modes  in  the  upstream  terrain.  This 
suggested  that  either  the  atmosphere  must  force  itself  into  one  of  these  wavelengths  or  the 
terrain  must  force  the  atmosphere  into  one  of  these  wavelengths  in  order  for  the  most 
severe  downslope  winds  to  occur  at  USAFA. 

Results  from  the  composite  studies  of  fourteen  severe  DSWS  link  the  atmospheric 
structures  found  in  NWB  events  to  the  structures  successfully  identified  by  ESA.  These 
results  show  that  ESA  successfully  identified  an  atmospheric  structure  that  is  common  in 
the  majority  of  NWB  DSWS  at  USAFA.  This  suggests  that  targeting  observations  at 
these  locations  prior  to  the  onset  of  a  non-breaking  severe  DSWS  would  impact  and 
potentially  improve  forecasts  of  these  phenomena  in  numerical  models. 
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C.  FUTURE  WORK 

Further  study  of  the  topics  put  forth  in  this  research  falls  into  two  main  areas: 
ESA  testing  and  the  further  investigation  into  the  initial  conditions  and  physical 
structures  that  lead  to  strong  DSWS,  flow  separation,  and  adjacent  lee  rotors  in  realistic 
complex  terrain. 

1.  ESA  Testing 

ESA  showed  skill  at  identifying  initial  conditions  that  affected  a  forecast  metric. 
However,  the  length  of  the  forecast  in  this  study  was  only  three  hours.  Future  studies 
should  focus  on  identifying  the  length  of  time  ESA  is  useful  for  highly  mesoscale  and 
complex  phenomena.  The  total  energy  perturbation  analysis  performed  on  the  ensemble 
runs  in  this  study  suggested  that  impacts  were  made  to  the  forecast  for  at  least  12  hours 
and  15-18  hours  in  isolated  cases.  The  length  of  the  impact  may  have  been  limited  in  this 
study  due  to  an  approaching  cold  front  that  changed  the  flow  regime.  Future  studies  may 
benefit  from  pushing  back  the  initialization  time  for  the  ensemble  runs  in  order  to  not  be 
as  close  to  the  end  of  the  downslope  event.  Doing  so  would  allow  for  the  testing  of  ESA 
with  longer  forecast  times.  Longer  forecast  times  also  may  allow  for  a  more  complete 
assessment  of  the  impact  of  ESA  on  ensemble  spread  over  time,  shedding  more  light  on 
whether  or  not  sensitivities  can  be  used  to  reduce  forecast  uncertainty  in  an  ensemble. 

ESA  was  applied  in  this  study  to  a  non- wave  breaking  severe  DSWS  and  found  to 
be  useful  tool.  Future  research  should  be  conducted  utilizing  ESA  and  WB  events  to 
determine  if  the  same  utility  applies.  Similar  methods  and  linearity  testing  should  be 
accomplished  to  address  this  unsettled  question. 

Additional  insight  utilizing  the  concepts  of  ESA  could  be  gained  by  switching  J 
and  x  in  ESA  and  doing  a  reverse  regression  (G.  Hakim,  personal  communication, 
December  12,  2012).  Instead  of  calculating  clJ/dx,  perhaps  calculating  dx/dJ  would  add 
insight  into  how  initial  conditions  (v)  impact  a  forecast  metric  ./.  This  would  pose  the 
question:  if  a  certain  J  is  wanted,  how  must  x  change  in  order  to  obtain  the  desired 
result?  An  experiment  of  this  type  could  allow  assessment  of  whether  one  sounding  at  a 
location  is  more  relevant  than  ten  surface  observations  in  a  different  location. 
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2.  Initial  Condition  Impacts  on  Downslope  Windstorms 

One  of  the  largest  impacts  of  this  study  is  the  assessment  that  the  leeside 
conditioning  plays  an  important  role  in  the  distance  DSWS  extend  down  the  slope,  the 
presence  and  location  of  flow  separation,  and  location  and  strength  of  adjacent  rotor’s 
return  flow.  Insight  into  the  role  of  leeside  conditioning  could  be  accomplished  using  a 
two-dimensional  flow  model  with  idealized  terrain.  This  would  allow  for  easy 
introduction  of  various  downstream  stability  and  wind  profiles,  and  a  fairly 
straightforward  method  for  evaluating  the  general  impact  of  the  leeside  environment.  In 
addition,  a  study  of  this  type  could  quantify  the  distance  downstream  perturbations  to  the 
leeside  environment  impact  upstream  downslope  winds,  flow  separation,  and  rotor 
formation.  Additionally,  a  study  of  this  type  could  help  determine  if  leeside  conditioning 
downstream  of  the  slope  plays  a  larger  role  in  determining  the  strength  of  a  DSWS  in 
non-wave  breaking  versus  wave  breaking  cases. 

Another  experiment  that  could  be  conducted  to  assess  the  impact  of  assimilating 
leeside  observations  would  be  to  assimilate  the  entire  or  smaller  subset  of  the  HWAS 
observations  into  an  ensemble  during  the  conditioning  period.  Ensemble  runs  at  various 
event  lead  times  and  lengths  could  be  run  to  quantify  the  impact  of  the  HWAS 
observations  on  the  forecast.  This  could  be  accomplished  in  conjunction  with  ESA 
predictions  to  access  ESA’s  ability  to  predict  the  location  of  initial  condition  sensitivity. 
A  study  of  this  type  would  help  quantify  a  mesonet’s  impact  on  modeling  and  forecast 
products  for  high  impact  events  in  complex  terrain. 

Additional  investigation  into  the  role  of  the  height  and  structure  of  the  tropopause 
on  DSWS,  flow  separation,  and  leeside  rotor  formation  is  needed.  Although  all 
perturbations  were  introduced  in  the  lower  troposphere  between  3  and  4  km  in  this  study, 
large  changes  in  the  temperature  and  wind  fields  were  seen  at  tropopause  level.  If  a 
change  in  the  lower  troposphere  caused  changes  at  the  tropopause  then  it  is  reasonable  to 
assume  that  sensitivities  exist  at  the  tropopause  level.  Previous  studies  have  noted  this 
sensitivity,  but  much  is  left  to  be  explored  in  this  area,  and  ESA  provides  a  good  tool  to 
assess  the  potential  impact  of  tropopause  height  and  structure  on  mountain  wave 
phenomena. 
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Lastly,  more  research  is  needed  to  enhance  our  understanding  of  the  role  of 
complex  terrain  upstream  of  the  DSWS  location.  Most  idealized  studies  in  the  past 
introduce  unperturbed  flow  upstream  of  regularly  shaped  terrain.  In  reality,  flow  is  rarely 
unperturbed  when  it  arrives  in  Colorado  due  to  the  multiple  upstream  mountain  ranges  in 
the  intermountain  west.  Studies  should  be  conducted  that  determine  the  wavelength  or 
frequency  of  perturbed  flow  before  it  encounters  the  Rocky  Mountains  and  how  those 
frequencies  relate  to  the  frequencies  of  highest  amplitude  found  by  the  Fourier  Analysis 
in  this  study.  It  may  be  found  that  perturbations  or  waves  already  inherent  in  the  upstream 
flow  play  an  important  role  in  determining  the  type  and  strength  of  DSWS  in  the  lee  of 
the  Rocky  Mountains  and  similar  locations  with  multiple  upstream  mountain  ranges. 
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APPENDIX  A 


A.  USAFA  GEOGRAPHY 

Appendix  A  follows  Homan  (2007). 

USAFA  lies  at  the  base  of  the  eastern  slopes  of  the  Front  Range  of  the  Rocky 
Mountains  approximately  20  km  northwest  of  Colorado  Springs,  Colorado.  The 
Continental  Divide  runs  north  to  south  through  Colorado  and  is  110  km  to  the  northwest 
of  USAFA  and  135  km  directly  to  the  west  (Figure  63).  The  average  ridgeline  elevation 
along  the  Continental  Divide  is  3600  m  with  several  peaks  greater  than  4200  m.  Although 
not  entirely  part  of  the  Continental  Divide,  the  Mosquito  Range  lies  about  105  km  to  the 
west-northwest  of  USAFA  with  ridgelines  of  just  over  3600  m  and  several  peaks  above 
4200  m.  The  Tarryall  and  Kenosha  Mountains  lie  even  closer  at  70-80  km  to  the 
northwest  with  several  peaks  in  the  3600  to  3800  m  range. 

The  local  dominant  feature  around  USAFA  is  the  Rampart  Range.  It  runs  north  to 
south  along  the  western  flank  of  USAFA  and  the  surrounding  areas.  Pikes  Peak  (4300  m) 
is  the  dominant  feature  of  the  Rampart  Range  and  is  located  approximately  21  km 
southwest  of  the  USAFA  Airfield.  The  average  ridgeline  of  the  Rampart  Range  directly 
adjacent  to  USAFA  is  about  2700  m  with  several  peaks  around  2900  m. 

The  USAFA  Reservation  is  18,500  acres  and  varies  greatly  in  elevation  both 
north  to  south  and  east  to  west.  There  are  several  ridges  that  extend  west  to  east  from  the 
Rampart  Range  along  USAFA,  with  their  average  elevations  ranging  from  2225  to  2000 
m.  Between  these  ridges  are  valleys  that  range  in  elevation  from  2050  to  1950  m.  The 
main  cadet  area  (often  referred  to  as  the  Terrazzo)  where  the  dormitories  and  academic 
buildings  are  located  is  positioned  on  one  of  these  ridges  at  2212  m.  The  USAFA  Airfield 
lies  at  the  mouth  of  Pine  and  Douglas  Valleys  with  elevations  ranging  from  1990  to  2005 
m.  Jack’s  Valley,  the  main  summertime  training  area  lies  on  the  north  end  of  the  USAFA 
Reservation  and  ranges  in  elevation  from  2130  m  to  the  west  and  1980  m  to  the  east. 
Monument  Creek  flows  north  to  south  along  the  eastern  edge  of  USAFA  and  parallels 
Interstate  25  (1-25)  and  a  set  of  railroad  tracks  from  the  northern  boundary  of  USAFA  to 
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roughly  where  1-25  is  adjacent  to  the  stadium.  At  this  point  the  creek  and  railroad  tracks 
turn  in  a  southwesterly  direction  and  are  west  of  the  USAFA  Airfield  where  the  creek 
turns  again  roughly  south  and  exits  the  USAFA  Reservation  just  southwest  of  the  South 
Gate.  Monument  Creek  is  the  center  of  a  valley  that  runs  north  to  south  with  elevations 
increasing  to  the  west  and  also  to  the  east. 


Figure  63.  Image  showing  mountain  ranges,  terrain  and  geography  of  Colorado.  Image 
taken  after  first  snowfall  in  2002.  (After  Descloitres,  2002). 

Located  about  18  km  north  of  the  USAFA  Airfield  and  13  km  north  of  the  cadet 
area  is  the  crest  of  the  Palmer/Monument  Divide  (sometimes  referred  to  as  the 
Palmer/Monument  Ridge).  This  ridge  extends  east  from  the  Rampart  Range  with  an 
average  elevation  of  2280  m  and  maximum  elevation  above  2400  m.  This  ridge  is  the 
main  geographical  feature  between  Colorado  Springs  and  Denver  and  measures  north  to 
south  on  average  45-55  km  and  runs  from  the  Rampart  Range  45-55  km  east  out  onto  the 
adjacent  plains.  Due  to  the  Palmer  Divide,  elevations  increase  as  you  head  north  from 
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USAFA  as  well  as  to  the  east.  The  region  directly  6-12  km  to  the  east  of  USAFA  on  the 
Palmer  Divide  is  known  as  the  Black  Forest  due  to  its  abundance  of  Ponderosa  Pines  with 
elevations  of  1980  m  to  2375  m.  The  terrain  drops  to  the  east  of  Black  Forest;  however, 
this  drop  is  much  more  pronounced  to  the  southeast  and  south  of  the  Academy.  The 
elevation  at  Limon,  CO  (87  km  to  the  east-northeast)  is  1635  m  while  at  Pueblo,  CO  (69 
km  to  the  south-southeast)  is  1440  m. 

B.  HIGH  WIND  ALERT  SYSTEM  (HWAS)  DESCRIPTION 

The  twelve  weather  observation  sensors  that  make  up  the  HWAS  network  are 
strategically  positioned  at  locations  to  provide  an  adequate  sample  of  meteorological  data  and 
give  a  good  representation  of  weather  conditions  around  USAFA  in  the  diverse  terrain.  A 
map  of  these  sensors  is  provided  in  Figure  64,  and  the  elevation  and  longitude/latitude 
positioning  of  each  sensor  is  provided  in  Table  8.  Nine  of  the  sensors  are  located  on  the 
USAFA  Reservation  with  the  Rampart  sensor  located  at  the  top  of  the  Rampart  Range 
ridgeline  due  west  of  USAFA  in  Pikes  Peak  National  Forest.  The  Pine  Creek  sensor  is 
located  on  a  roof  at  a  local  school  east  of  USAFA  while  the  Lewis  Palmer  sensor  (also  on  a 
roof  at  a  local  school)  is  north  of  USAFA.  Both  of  these  sensors  are  on  the  edge  of  the 
Palmer  Divide  and  the  terrain  decreases  in  elevation  from  these  sensors  toward  USAFA. 


Sensor  Name 

Latitude 

Longitude 

Elevation  (m) 

Aardvark 

39°  02’  07”  N 

104°  50’  39” 

W 

2042.2 

Airfield 

38°  56’  06”  N 

104°  49’  30” 

W 

1978.2 

Command  Post 

39°  00’  37”  N 

104°  53’  11” 

W 

2179.3 

Community 

Center 

38°  58’  57”  N 

104°  52’  29” 

W 

2185.4 

Lewis  Palmer 

39°  04’  48”  N 

104°  51’  16” 

W 

2132.1 

North  Ridge 

39°  01’  33”  N 

104°  52’  43” 

W 

2183.9 

Pine  Creek 

38°  59’  12”  N 

104°  45’  48” 

W 

2126.0 

Pine  Valley 

38°  58’  13”  N 

104°  51’  37” 

W 

2026.9 

Rampart 

38°  59’  13”  N 

104°  55’  14” 

W 

2840.8 

South  Gate 

38°  56’  56”  N 

104°  48’  55” 

W 

1953.8 

South  Ridge 

38°  57’  25”  N 

104°50’  51” 

W 

2036.1 

Stadium 

38°  59’  51”  N 

104°  50’  55” 

W 

2084.9 

Table  8.  HWAS  sensors  longitude/latitude  and  elevation. 
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HWAS  sensors  located  on  the  USAFA  Reservation  are  mostly  located  along  west 
to  east  oriented  ridgelines  with  a  few  exceptions.  The  Pine  Valley  sensor  is  on  the  rooftop 
of  a  school  that  is  located  in  the  largest  west  to  east  oriented  valley  on  USAFA.  The 
Airfield  sensor  is  on  the  east  side  of  Monument  Creek  and  is  located  at  the  mouth  of  Pine 
Valley.  This  sensor  is  in  an  open  field  just  west  of  the  adjacent  runways  on  the  USAFA 
Airfield  for  the  case  studies  in  this  research,  but  was  moved  to  a  location  northwest  of  the 
mnways  in  2011.  The  Aardvark  sensor  lies  along  the  eastern  most  extension  of  the 
Northridge.  It  is  adjacent  to  an  additional  runway  used  for  training  and  Monument  Creek. 
When  strong  winds  are  present  from  the  west,  because  this  sensor  is  not  in  a  west-east 
oriented  valley  but  at  the  base  of  a  ridgeline,  this  sensor  is  often  spared  the  strongest 
winds  and  the  recorded  wind  values  are  much  lower  than  at  other  sensors.  Instead,  with 
its  proximity  to  Monument  Creek,  it  is  in  an  ideal  location  to  measure  drainage  flow  from 
the  north  coming  down  the  Palmer  Divide.  The  only  other  sensor  that  is  located  on  a 
platform  that  is  on  top  of  a  building  is  the  Command  Post  sensor  which  is  located  on  top 
of  Vandenberg  Hall,  a  six-story  tall  cadet  dormitory.  This  sensor  has  the  greatest  number 
of  man-made  obstructions  around  it  of  any  of  the  sensors  with  several  other  buildings 
close  by  in  all  directions.  It  is  the  furthest  west  and  closest  to  the  steep  slopes  of  the 
Rampart  Range  than  any  of  the  other  sensors  except  for  Rampart  which  lies  at  the  atop  of 
the  Front  Range. 
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Figure  64.  Map  of  locations  of  HWAS  sensors.  Green  shaded  regions  indicate  areas 

covered  by  forest. 
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APPENDIX  B 


A.  PERFECT  MODEL— “TRUTH  RUN”  FIGURES 
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Figure  65.  Terrain  height  in  meters  for  innermost  domain  centered  over  Colorado.  Solid 
black  hne  indicates  location  of  cross  sections  shown  in  upper  plots  of  figures 
that  follow  at  the  gridline  of  the  Rampart  sensor. 
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Figure  66.  Top:  Cross  section  of  x-wind  component  (shaded)  and  potential  temperature 
(contoured)  from  west  to  east  along  solid  black  line  in  previous  figure  for  the 
12/30/2008  00Z  truth  run.  Black  inset  box  indicates  location  of  zoomed  in 
cross  section  depicted  at  bottom  half  of  figure.  The  horizontal  axis  is  scaled 
to  the  number  of  grid  points  from  the  innermost  domain  with  a  resolution  of 
1.33  km.  Approximate  locations  of  four  HWAS  sensors  are  marked  by 

circles. 
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Figure  67.  Truth  run  for  12/30/2008  03Z.  Description  follows  previous  figure. 
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Figure  68.  Truth  run  for  12/30/2008  06Z.  Description  follows  previous  figure. 
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Figure  69.  Truth  run  for  12/30/2008  09Z.  Description  follows  previous  figure. 
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Figure  70.  Truth  run  for  12/30/2008  12Z.  Description  follows  previous  figure. 
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Figure  71.  Truth  run  for  12/30/2008  15Z.  Description  follows  previous  figure. 
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Figure  72.  Truth  run  for  12/30/2008  18Z.  Description  follows  previous  figure. 
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Figure  73.  Truth  run  for  12/30/2008  21Z.  Description  follows  previous  figure. 
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Figure  74.  Truth  run  for  12/31/2008  00Z.  Description  follows  previous  figure. 
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APPENDIX  C 


A.  FIGURES  OF  WRF  MODEF  RUNS  FOR  FOURTEEN  SEVERE 

DOWNSLOPE  WINDSTORMS  AT  USAFA 

Appendix  C  contains  cross  sections  for  each  of  the  fourteen  severe  DSWS  for  the 
case-study  presented  in  Chapter  V.  The  time  of  each  cross  section  corresponds  to  the  time 
of  maximum  winds  as  predicted  by  the  model  in  the  solid  black  box  over  USAFA  shown 
in  the  top  of  each  figure  and  shown  in  the  zoomed-in  cross  section  in  the  bottom  of  each 
figure.  The  figures  are  ordered  with  the  six  non-breaking  events  shown  first  and  then  the 
eight  vertically  propagating  breaking  wave  events  shown  second.  The  times  displayed  in 
these  figures  are  not  necessarily  the  best  example  of  larger  scale  trapping  or  wave¬ 
breaking  used  to  classify  each  event,  but  rather  correspond  to  the  time  of  maximum 
winds  lee  of  the  Rampart  Range  at  USAFA. 
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Figure  75.  Top:  Cross  section  from  west  to  east  as  depicted  by  solid  black  line  in  Figure 
16  for  the  1 1/03/2005  17Z  WRF  run.  Black  inset  box  indicates  location  of 
zoomed  in  cross  section  depicted  at  bottom  half  of  figure.  The  x-wind 
component  is  shaded  in  color  and  theta  is  contoured  with  solid  black  lines. 
The  horizontal  axis  is  scaled  to  the  number  of  grid  points  from  the  innermost 
domain  with  a  resolution  of  1.33  km.  Approximate  locations  of  four  HWAS 

sensors  are  marked  by  circles. 
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Figure  76.  WRF  run  for  03/06/2004  18Z.  Description  follows  previous  figure. 
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Figure  77.  WRF  truth  run  for  12/30/2008  12Z.  Description  follows  previous  figure. 
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Figure  78.  WRF  run  for  01/19/2012  04Z.  Description  follows  previous  figure. 
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Figure  79.  WRF  run  for  12/20/2010  05Z.  Description  follows  previous  figure. 
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Figure  80.  WRF  run  for  11/12/2011  18Z.  Description  follows  previous  figure. 
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Figure  81.  WRF  run  for  12/18/2008  23Z.  Description  follows  previous  figure. 
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Figure  82.  WRF  run  for  01/05/2008  15Z.  Description  follows  previous  figure. 
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Figure  83.  WRF  run  for  02/21/2012  12Z.  Description  follows  previous  figure. 
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Figure  85.  WRF  run  for  10/02/2007  09Z.  Description  follows  previous  figure. 
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Figure  86.  WRF  run  for  1 1/13/201 1  07Z.  Description  follows  previous  figure. 
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Figure  87.  WRF  run  for  06/07/2007  13Z.  Description  follows  previous  figure. 
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Figure  88.  WRF  run  for  04/15/2006  13Z.  Description  follows  previous  figure. 


149 


THIS  PAGE  INTENTIONALLY  LEFT  BLANK 


150 


LIST  OF  REFERENCES 


Aanenson,  C.  J.  M.,  1965:  Gales  in  Yorkshire  in  February  1962.  Geophys.  Mem.,  14, 
1-44. 

Ancell,  B.,  &  Hakim,  G.,  2007:  Comparing  adjoint-  and  ensemble-sensitivity  analysis 
with  applications  to  observation  targeting.  Mon.  Wea.  Rev.,  135,  4117-4134. 

Anderson,  J.  L.,  2003:  A  local  least  squares  framework  for  ensemble  filtering.  Mon.  Wea. 
Rev.,  131,  634-642 

Anderson,  J.,  T.  Hoar,  K.  Raeder,  H.  Liu,  N.  Collins,  R.  Torn,  and  A.  Avellano,  2009: 

The  data  assimilation  research  testbed:  A  community  facility.  Bull.  Amer. 

Meteor.  Soc.,  90,  1283-1296. 

Berliner,  L.  M.,  Lu,  Z.-Q.,  &  Snyder,  C.,  1999:  Statistical  design  for  adaptive  weather 
observations.  J.  Atmos.  Sci.,  56,  2536-2552. 

Brinkman,  W.,  1974:  Strong  downslope  winds  at  boulder,  Colorado.  Mon.  Wea.  Rev., 

102,  592-602. 

Clark,  T.  L.,  and  R.  D.  Farley,  1984:  Severe  DSWS  calculations  in  two  and  three  spatial 
dimensions  using  an  elastic  interactive  grid  nesting:  A  possible  mechanism  for 
gustiness.  J.  Atmos.  Sci.,  41,  329-350. 

Clark,  T.  L.,  and  W.  R.  Peltier,  1977:  On  the  evolution  and  stability  of  finite  amplitude 
mountain  waves.  J.  Atmos.  Sci.,  34,  1715-1730. 

Clark,  T.  L.  and  Peltier,  W.  R.,  1984:  Critical  level  reflection  and  resonant  growth  of 
nonlinear  mountain  waves.  J.  Atmos.  Sci.,  41,  3122-3134. 

Descloitres,  J.,  MODIS  Land  Rapid  Response  Team,  cited  March  2007.  [Available 
online  at  http://www.14ers.eom/images/ranges3db700x414.ipg.l 

Doyle,  J.,  Amerault,  C.,  &  Reynolds,  C.  A.,  2007:  Sensitivity  analysis  of  mountain 
waves  using  an  adjoint  model.  Meteorologische  Zeitschrift,  16,  607-620. 

Doyle,  J.  D.,  and  D.  R.  Durran,  2002:  The  dynamics  of  mountain-wave-induced  rotors.  J. 
Atmos.  Sci.,  59,  186-201. 

Durran,  D.R.,  1986:  Another  look  at  downslope  winds.  Part  I:  The  development  of 

analogs  to  supercritical  flow  in  an  infinitely  deep,  continuously  stratified  fluid.  J. 
Atmos.  Sci.,  43,  2527-2543. 

Durran,  D.  R.,  1990:  Mountain  waves  and  downslope  winds.  Atmospheric  Processes 
Over  Complex  Terrain.  Amer.  Meteor.  Soc.,  23,  59-81. 

151 


Durran,  D.R.,  2003:  Downslope  winds.  Encyclopedia  of  atmospheric  sciences,  2003, 
Holton,  J.R.,  J.  Pyle  and  J.A.  Curry  eds.  pp.  1161-1170,  Elsevier  Science  Ltd. 

Durran,  D.  R.,  and  J.  B.  Klemp,  1987:  Another  look  at  downslope  winds.  Part  II: 

Nonlinear  amplification  beneath  wave-overturning  layers.  J.  Atmos.  Sci.,  22, 
3402-3412. 

Earth  System  Research  Laboratory,  cited  July  2013:  Boulder  wind  info. 

[Available  online  at  ESRL  :  PSD  :  Boulder  Wind  Events.l 

Eliassen,  A.,  and  E.  Palm,  1960:  On  the  transfer  of  energy  in  stationary  mountain  wave. 
Geofis.  Pub.,  22,  1-23. 

Gaspari,  G.,  and  S.  E.  Cohn,  1999:  Construction  of  correlation  functions  in  two  and  three 
dimensions.  Quarterly  Journal  of  the  Royal  Meteorological  Society,  125,  723- 
757. 

Glickman,  T.  S.,  Ed.,  2000:  Glossary  of  Meteorology.  2d  ed.  Amer.  Meteor.  Soc.,  855  pp. 

Grimit,  E.  P.,  and  C.  F.  Mass,  2007:  Measuring  the  ensemble  spread-error  relationship 
with  a  probabilistic  approach:  Stochastic  ensemble  results.  Mon.  Wea.  Rev.,  135, 
203-221. 

Hakim,  G.,  &  Torn,  R.,  2008:  Ensemble  synoptic  analysis.  Meteorological  Monographs, 
33,  147-161. 

Hamill,  T.  M,  cited  2013:  Review  of  “covariance  localization”  in  ensemble  filters. 
[Available  online  at 

http://www.esrl.noaa.gov/psd/people/tom.hamill/CovLocal  BuenosAires  Hamill. 

ppt.] 

Hertenstein,  R.  F.,  2009:  The  influence  of  inversions  on  rotors.  Mon.  Wea.  Rev.,  137, 
433-446. 

Hertenstein,  R.  F.,  and  J.  P.  Kuettner,  2005:  Rotor  types  associated  with  steep  lee 
topography:  Influence  of  the  wind  profile.  Tellus,  51  A,  117-135. 

Homan,  P.,  2007:  Evaluation  of  high  density  surface  observations  in  complex  terrain  and 
their  contribution  to  the  MM5  model.  Naval  Postgraduate  School. 

Holmboe,  J.,  and  H.  Klieforth,  1957:  Investigations  of  mountain  lee  waves  and  airflow 
over  the  Sierra  Nevada.  Final  Rep.  Contract  AF19(604)-728,  University  of 
California  ADNo.  133606,  Department  of  Meteorology,  University  of  California, 
Los  Angeles,  290  pp. 

Houghton,  D.  D.,  and  E.  Isaacson,  1968:  Mountain  winds.  Stud.  Num.  Anal.,  2,  21-52. 


152 


Houtekamer,  P.  L.,  and  Mitchell,  H.  L,  2001:  A  sequential  Kalman  Ensemble  Filter  for 
atmospheric  data  assimilation.  Mon.  Wea.  Rev.,  129,  123-137. 

Jiang,  Q.,  and  Doyle,  J.,  2008:  Diurnal  variation  of  downslope  winds  in  Owens  Valley 
during  the  Sierra  rotor  experiment.  Mon.  Wea.  Rev..  136,  3760-3780. 

Jiang, Q.,  J.  D.  Doyle,  S.  Wang,  and  R.  B.  Smith,  2007:  On  boundary  layer  separation  in 
the  lee  of  mesoscale  topography.  J.  Atmos.  Sci.,  64,  401^420. 

Kalnay,  E.,  2009:  Atmospheric  Modeling,  Data  Assimilation  and  Predictability .  5th  ed. 
Cambridge  University  Press,  341  pp. 

Kim,  H.  M.,  and  Jung,  B.  J.,  2009:  Influence  of  moist  physics  and  norms  on  singular 
vectors  for  a  tropical  cyclone.  Mon.  Wea.  Rev.,  137,  525-543. 

Klemp,  J.B.,  and  Lilly,  D.R.,  1975:  The  dynamics  of  wave-induced  downslope  winds.  J. 
Atmos.  Sci.,  32,  320-339. 

Kuettner,  J.,  1959:  The  rotor  flow  in  the  lee  of  mountains.  GRD  Res.  Notes  No.  6, 
Geophys.  Res.  Directorate,  Air  Force  Cambridge  Res.  Center. 

Lee,  T.,  Rielke,  R.,  Kessker,  R.,  &  Weaver,  J.,  1989:  Influence  of  cold  pools  downstream 
of  mountain  barriers  on  downslope  winds  and  flushing.  Mon.  Wea.  Rev.,  117, 
2041-2058. 

Lester,  P.  F.,  and  Fingergut,  1974:  Lower  turbulent  zones  associated  with  mountain  lee 
waves.  J.  Appl.  Meteor.,  13,  54-61. 

Lilliefors,  H.  W.,  1967:  On  the  Kolmogorov-Smirnov  test  for  normality  with  mean  and 
variance  unknown.  Journal  of  the  American  Statistical  Association.  62, 

399-402. 

Lorenz,  Edward  N.,  1963:  Deterministic  nonperiodic  flow.  J.  Atmos.  Sci.,  20,  130-141. 

Lorenz,  E.  N.,  1969:  The  predictability  of  a  flow  which  contains  many  scales  of  motion. 
Tellus,  21,  289-307. 

Lorenz,  E.,  Emanuel,  K.,  1998:  Optimal  sites  for  supplementary  weather  observations: 
simulation  with  a  small  model.  J.  Atmos.  Sci.,  55,  399-414. 

Long,  R.  R.,  1953a:  Some  aspects  of  the  flow  of  stratified  fluids.  I.  A  theoretical 
investigation.  Tellus,  5,  42-58. 

METED/COMET  Program,  cited  July  2013:  Mountain  Waves  and  Downslope  Winds 
rhttp://www.meted.ucar.edu/mesoprim/mtnwave/print.htm.l 


153 


Oard,  Michael  J.,  1993:  A  method  for  predicting  Chinook  winds  east  of  the  Montana 
Rockies.  Wea.  Forecasting,  8,  166-180. 

Peltier,  W.  R.  and  T.  L.  Clark,  1979:  The  evolution  and  stability  of  finite  amplitude 
mountain  waves  II:  Surface  wave  drag  and  severe  DSWS.  J.  Atmos.  Sci.,  36, 
1498-1529. 

Pielke,  R.A.,  1985:  The  use  of  mesoscale  numerical  models  to  assess  wind  distribution 
and  boundary  layer  structure  in  complex  terrain.  Bound. -Layer  Meteor.,  31, 
217-231. 

Reinecke,  P.,  and  Durran,  D.,  2008:  Estimating  topographic  blocking  using  a  Froude 
number  when  the  static  stability  is  non-uniform.  J.  Atmos.  Sci.,  65,  1035-1048 

Reinecke,  P.,  and  Durran,  D.,  2009:  Initial-condition  sensitivities  and  the  predictability 
of  downslope  winds.  J.  Atmos.  Sci.,  66,  3401-3418. 

Rocky  Mountain  Insurance  Information  Association,  cited  January  2009:  Tuesday’s 

gusts  blow  in  Colorado’s  fourth  most  expensive  windstorm.  [Available  online  at 
http://www.rmiia.org/News  room/catastrophe  news/2009  01  02  Colorado  Win 

d.asp.l 

Saha,  et  al.,  2010:  The  NCEP  Climate  Forecast  System  Reanalysis.  Bull.  Amer.  Meteor. 
Soc.,  91,  1015-1057. 

Skamarock,  W.C.,  J.  B.  Klemp,  J.  Dudhia,  D.  O.  Gill,  D.  M.  Barker,  M.  G.  Duda,  X. 
Huang,  W.  Wang,  and  J.  G.  Powers,  2008:  A  description  of  the  advanced 
research  WRF  version  3.  NCAR  Tech.  Note,  NCAR/TN-475+STR,  113  pp. 

Scorer,  R.S.,  and  H.  Klieforth,  1959:  Theory  of  mountain  waves  of  large  amplitude. 
Quart.  J.  Roy.  Meteor.  Soc.,  85,  131-143. 

Smith,  R.  B.,  1977:  The  steepening  of  hydrostatic  mountain  waves.  J.  Atmos.  Sci.,  34, 
1634-1654. 

Snyder,  C.,  1996:  Summary  of  an  informal  workshop  on  adaptive  observations  and 
FASTEX.  Bull.  Amer.  Meteor.  Soc.,  77,  953-961. 

Stull,  R.  B.,  1988:  An  Introduction  to  Boundary  Layer  Meteorology.  Kluwer  Academic 
Publishers,  666. 

Torn,  R.,  2010:  Ensemble -based  sensitivity  analysis  applied  to  African  easterly  waves. 
Weather  and  Forecasting,  25,  61-78. 

Torn,  R.,  &  Hakim,  G.,  2008:  Ensemble-based  sensitivity  analysis.  Mon.  Wea.  Rev.,  136, 
663-677. 


154 


Torn,  R.,  and  Hakim,  G.,  2009:  Initial  condition  sensitivity  of  Western  Pacific 

extratropical  transitions  determined  using  ensemble-based  sensitivity  analysis. 
Mon.  Wea.  Rev.,  137,  3388-3406. 

Torn,  R.,  G.  Hakim  and  C.  Snyder,  2006:  Boundary  conditions  for  limited-area  ensemble 
Kalman  filters.  Mon.  Wea.  Rev.,  134,  2490-2502. 

The  Weather  Research  and  Forecasting  Model  (WRF),  cited  2013:  WRF  Model  Users 
Page.  [Available  online  at  http://wrf-model.Org/users/users.php.l 

Vosper,  S.  B.,  2004:  Inversion  effects  on  mountain  lee  waves.  Quart.  J.  Roy.  Meteor. 
Soc.,  130,  1723-1748. 

Vosper,  S.  G.,  P.  F.  Sheridan,  and  A.  R.  Brown,  2006:  Flow  separation  and  rotor 

formation  beneath  two-dimensional  trapped  lee  waves.  Quart.  J.  Roy.  Meteor. 

Soc.,  132,  2415-2438. 

Wilks,  D.  S.,  2006:  Statistical  methods  in  the  atmospheric  sciences.  2nd  ed.  Academic 
Press,  627  pp. 

Young,  G.,  and  Pielke,  R.,  1983:  Application  of  terrain  height  variance  spectra  to 
mesoscale  modeling.  J.  Atmos.  Sci.,  40,  2555-2560. 

Zou,  X.,  F.  Vandenberghe,  M.  Pondeca,  and  Y.-H.  Kuo,  1997:  Introduction  to  adjoint 

techniques  and  the  MM5  adjoint  modeling  system.  NCAR  Tech.  Note  NCAR/TN- 
435STR,  110  pp. 


155 


THIS  PAGE  INTENTIONALLY  LEFT  BLANK 


156 


INITIAL  DISTRIBUTION  LIST 


1.  Defense  Technical  Information  Center 
Ft.  Belvoir,  Virginia 

2.  Dudley  Knox  Library 
Naval  Postgraduate  School 
Monterey,  California 


157 


